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

Diff for /OpenXM_contrib2/asir2000/engine/NEZ.c between version 1.1 and 1.7

version 1.1, 1999/12/03 07:39:08 version 1.7, 2018/03/29 01:32:51
Line 1 
Line 1 
 /* $OpenXM: OpenXM/src/asir99/engine/NEZ.c,v 1.1.1.1 1999/11/10 08:12:26 noro Exp $ */  /*
    * Copyright (c) 1994-2000 FUJITSU LABORATORIES LIMITED
    * All rights reserved.
    *
    * FUJITSU LABORATORIES LIMITED ("FLL") hereby grants you a limited,
    * non-exclusive and royalty-free license to use, copy, modify and
    * redistribute, solely for non-commercial and non-profit purposes, the
    * computer program, "Risa/Asir" ("SOFTWARE"), subject to the terms and
    * conditions of this Agreement. For the avoidance of doubt, you acquire
    * only a limited right to use the SOFTWARE hereunder, and FLL or any
    * third party developer retains all rights, including but not limited to
    * copyrights, in and to the SOFTWARE.
    *
    * (1) FLL does not grant you a license in any way for commercial
    * purposes. You may use the SOFTWARE only for non-commercial and
    * non-profit purposes only, such as academic, research and internal
    * business use.
    * (2) The SOFTWARE is protected by the Copyright Law of Japan and
    * international copyright treaties. If you make copies of the SOFTWARE,
    * with or without modification, as permitted hereunder, you shall affix
    * to all such copies of the SOFTWARE the above copyright notice.
    * (3) An explicit reference to this SOFTWARE and its copyright owner
    * shall be made on your publication or presentation in any form of the
    * results obtained by use of the SOFTWARE.
    * (4) In the event that you modify the SOFTWARE, you shall notify FLL by
    * e-mail at risa-admin@sec.flab.fujitsu.co.jp of the detailed specification
    * for such modification or the source code of the modified part of the
    * SOFTWARE.
    *
    * THE SOFTWARE IS PROVIDED AS IS WITHOUT ANY WARRANTY OF ANY KIND. FLL
    * MAKES ABSOLUTELY NO WARRANTIES, EXPRESSED, IMPLIED OR STATUTORY, AND
    * EXPRESSLY DISCLAIMS ANY IMPLIED WARRANTY OF MERCHANTABILITY, FITNESS
    * FOR A PARTICULAR PURPOSE OR NONINFRINGEMENT OF THIRD PARTIES'
    * RIGHTS. NO FLL DEALER, AGENT, EMPLOYEES IS AUTHORIZED TO MAKE ANY
    * MODIFICATIONS, EXTENSIONS, OR ADDITIONS TO THIS WARRANTY.
    * UNDER NO CIRCUMSTANCES AND UNDER NO LEGAL THEORY, TORT, CONTRACT,
    * OR OTHERWISE, SHALL FLL BE LIABLE TO YOU OR ANY OTHER PERSON FOR ANY
    * DIRECT, INDIRECT, SPECIAL, INCIDENTAL, PUNITIVE OR CONSEQUENTIAL
    * DAMAGES OF ANY CHARACTER, INCLUDING, WITHOUT LIMITATION, DAMAGES
    * ARISING OUT OF OR RELATING TO THE SOFTWARE OR THIS AGREEMENT, DAMAGES
    * FOR LOSS OF GOODWILL, WORK STOPPAGE, OR LOSS OF DATA, OR FOR ANY
    * DAMAGES, EVEN IF FLL SHALL HAVE BEEN INFORMED OF THE POSSIBILITY OF
    * SUCH DAMAGES, OR FOR ANY CLAIM BY ANY OTHER PARTY. EVEN IF A PART
    * OF THE SOFTWARE HAS BEEN DEVELOPED BY A THIRD PARTY, THE THIRD PARTY
    * DEVELOPER SHALL HAVE NO LIABILITY IN CONNECTION WITH THE USE,
    * PERFORMANCE OR NON-PERFORMANCE OF THE SOFTWARE.
    *
    * $OpenXM: OpenXM_contrib2/asir2000/engine/NEZ.c,v 1.6 2010/02/19 08:27:19 noro Exp $
   */
 #include "ca.h"  #include "ca.h"
   
 void nezgcdnpz(vl,ps,m,pr)  void nezgcdnpz(vl,ps,m,pr)
Line 6  VL vl;
Line 54  VL vl;
 P *ps,*pr;  P *ps,*pr;
 int m;  int m;
 {  {
         P t,s,mg;    P t,s,mg;
         VL tvl,svl,avl,nvl;    VL tvl,svl,avl,nvl;
         int i,j,k;    int i,j,k;
         N c;    N c;
         Q cq;    Q cq;
         P *pl,*pm;    P *pl,*pm;
         DCP *ml;    DCP *ml;
         Q *cl;    Q *cl;
         P **cp;    P **cp;
         int *cn;    int *cn;
   
         pl = (P *)ALLOCA(m*sizeof(P));    if ( m == 1 ) {
         ml = (DCP *)ALLOCA(m*sizeof(DCP));      *pr = ps[0]; return;
         cl = (Q *)ALLOCA(m*sizeof(Q));    }
         for ( i = 0; i < m; i++ )    pl = (P *)ALLOCA(m*sizeof(P));
                 monomialfctr(vl,ps[i],&pl[i],&ml[i]);    ml = (DCP *)ALLOCA(m*sizeof(DCP));
         gcdmonomial(vl,ml,m,&mg);    cl = (Q *)ALLOCA(m*sizeof(Q));
         for ( i = 0; i < m; i++ ) {    for ( i = 0; i < m; i++ )
                 ptozp(pl[i],1,&cl[i],&t); pl[i] = t;      monomialfctr(vl,ps[i],&pl[i],&ml[i]);
         }    gcdmonomial(vl,ml,m,&mg);
         for ( i = 1, cq = cl[0]; i < m; i++ ) {    for ( i = 0; i < m; i++ ) {
                 gcdn(NM(cl[i]),NM(cq),&c); NTOQ(c,1,cq);      ptozp(pl[i],1,&cl[i],&t); pl[i] = t;
         }    }
         for ( i = 0; i < m; i++ )    for ( i = 1, cq = cl[0]; i < m; i++ ) {
                 if ( NUM(pl[i]) ) {      gcdn(NM(cl[i]),NM(cq),&c); NTOQ(c,1,cq);
                         mulp(vl,(P)cq,mg,pr); return;    }
                 }    for ( i = 0; i < m; i++ )
         for ( i = 0, nvl = vl, avl = 0; nvl && i < m; i++ ) {      if ( NUM(pl[i]) ) {
                 clctv(vl,pl[i],&tvl);        mulp(vl,(P)cq,mg,pr); return;
                 intersectv(nvl,tvl,&svl); nvl = svl;      }
                 mergev(vl,avl,tvl,&svl); avl = svl;    for ( i = 0, nvl = vl, avl = 0; nvl && i < m; i++ ) {
         }      clctv(vl,pl[i],&tvl);
         if ( !nvl ) {      intersectv(nvl,tvl,&svl); nvl = svl;
                 mulp(vl,(P)cq,mg,pr); return;      mergev(vl,avl,tvl,&svl); avl = svl;
         }    }
         if ( !NEXT(avl) ) {    if ( !nvl ) {
                 nuezgcdnpzmain(vl,pl,m,&t); mulp(vl,t,(P)cq,&s); mulp(vl,s,mg,pr);      mulp(vl,(P)cq,mg,pr); return;
                 return;    }
         }    if ( !NEXT(avl) ) {
         for ( tvl = nvl, i = 0; tvl; tvl = NEXT(tvl), i++ );      nuezgcdnpzmain(vl,pl,m,&t); mulp(vl,t,(P)cq,&s); mulp(vl,s,mg,pr);
         for ( tvl = avl, j = 0; tvl; tvl = NEXT(tvl), j++ );      return;
         if ( i == j ) {    }
                 /* all the pl[i]'s have the same variables */    for ( tvl = nvl, i = 0; tvl; tvl = NEXT(tvl), i++ );
                 sortplistbyhomdeg(pl,m); nezgcdnpzmain(nvl,pl,m,&t);    for ( tvl = avl, j = 0; tvl; tvl = NEXT(tvl), j++ );
     if ( i == j ) {
       /* all the pl[i]'s have the same variables */
       sortplistbyhomdeg(pl,m); nezgcdnpzmain(nvl,pl,m,&t);
 #if 0  #if 0
                 /* search the minimal degree poly */      /* search the minimal degree poly */
                 for ( i = 0; i < m; i++ ) {      for ( i = 0; i < m; i++ ) {
                 for ( tvl = nvl; tvl; tvl = NEXT(tvl) ) {            for ( tvl = nvl; tvl; tvl = NEXT(tvl) ) {
                 dt = getdeg(tvl->v,pl[i]);                  dt = getdeg(tvl->v,pl[i]);
                 if ( tvl == nvl || dt < d1 ) {                  if ( tvl == nvl || dt < d1 ) {
                         v1 = tvl->v; d1 = dt;                          v1 = tvl->v; d1 = dt;
                 }                  }
                 }            }
                         if ( i == 0 || d1 < d ) {        if ( i == 0 || d1 < d ) {
                                 v = v1; d = d1; j = i;          v = v1; d = d1; j = i;
                         }        }
                 }      }
                 t = pl[0]; pl[0] = pl[j]; pl[j] = t;      t = pl[0]; pl[0] = pl[j]; pl[j] = t;
                 if ( v != nvl->v ) {      if ( v != nvl->v ) {
                         reordvar(nvl,v,&mvl);        reordvar(nvl,v,&mvl);
                         for ( i = 0; i < m; i++ ) {        for ( i = 0; i < m; i++ ) {
                                 reorderp(mvl,nvl,pl[i],&t); pl[i] = t;          reorderp(mvl,nvl,pl[i],&t); pl[i] = t;
                         }        }
                         nezgcdnpzmain(mvl,pl,m,&s); reorderp(nvl,mvl,s,&t);        nezgcdnpzmain(mvl,pl,m,&s); reorderp(nvl,mvl,s,&t);
                 } else      } else
                         nezgcdnpzmain(nvl,pl,m,&t);        nezgcdnpzmain(nvl,pl,m,&t);
 #endif  #endif
         } else {    } else {
                 cp = (P **)ALLOCA(m*sizeof(P *));      cp = (P **)ALLOCA(m*sizeof(P *));
                 cn = (int *)ALLOCA(m*sizeof(int));      cn = (int *)ALLOCA(m*sizeof(int));
                 for ( i = 0; i < m; i++ ) {      for ( i = 0; i < m; i++ ) {
                         cp[i] = (P *)ALLOCA(lengthp(pl[i])*sizeof(P));        cp[i] = (P *)ALLOCA(lengthp(pl[i])*sizeof(P));
                         cn[i] = pcoef(vl,nvl,pl[i],cp[i]);        cn[i] = pcoef(vl,nvl,pl[i],cp[i]);
                 }      }
                 for ( i = j = 0; i < m; i++ )      for ( i = j = 0; i < m; i++ )
                         j += cn[i];        j += cn[i];
                 pm = (P *)ALLOCA(j*sizeof(P));      pm = (P *)ALLOCA(j*sizeof(P));
                 for ( i = k = 0; i < m; i++ )      for ( i = k = 0; i < m; i++ )
                         for ( j = 0; j < cn[i]; j++ )        for ( j = 0; j < cn[i]; j++ )
                                 pm[k++] = cp[i][j];          pm[k++] = cp[i][j];
                 nezgcdnpz(vl,pm,k,&t);      nezgcdnpz(vl,pm,k,&t);
         }    }
         mulp(vl,t,(P)cq,&s); mulp(vl,s,mg,pr);    mulp(vl,t,(P)cq,&s); mulp(vl,s,mg,pr);
 }  }
   
 void sortplistbyhomdeg(p,n)  void sortplistbyhomdeg(p,n)
 P *p;  P *p;
 int n;  int n;
 {  {
         int i,j,k;    int i,j,k;
         int *l;    int *l;
         P t;    P t;
   
         l = (int *)ALLOCA(n*sizeof(int));    l = (int *)ALLOCA(n*sizeof(int));
         for ( i = 0; i < n; i++ )    for ( i = 0; i < n; i++ )
                 l[i] = homdeg(p[i]);      l[i] = homdeg(p[i]);
         for ( i = 0; i < n; i++ )    for ( i = 0; i < n; i++ )
                 for ( j = i + 1; j < n; j++ )      for ( j = i + 1; j < n; j++ )
                         if ( l[j] < l[i] ) {        if ( l[j] < l[i] ) {
                                 t = p[i]; p[i] = p[j]; p[j] = t;          t = p[i]; p[i] = p[j]; p[j] = t;
                                 k = l[i]; l[i] = l[j]; l[j] = k;          k = l[i]; l[i] = l[j]; l[j] = k;
                         }        }
 }  }
   
 void nuezgcdnpzmain(vl,ps,m,r)  void nuezgcdnpzmain(vl,ps,m,r)
Line 116  P *ps;
Line 167  P *ps;
 int m;  int m;
 P *r;  P *r;
 {  {
         P *tps;    P *tps;
         P f,t;    P f,t;
         int i;    int i;
   
         for ( i = 0; i < m; i++ )    for ( i = 0; i < m; i++ )
                 if ( NUM(ps[i]) ) {      if ( NUM(ps[i]) ) {
                         *r = (P)ONE; return;        *r = (P)ONE; return;
                 }      }
         tps = (P *) ALLOCA(m*sizeof(P));    tps = (P *) ALLOCA(m*sizeof(P));
         for ( i = 0; i < m; i++ )    for ( i = 0; i < m; i++ )
                 tps[i] = ps[i];      tps[i] = ps[i];
         sortplist(tps,m);    sortplist(tps,m);
         for ( i = 1, f = tps[0]; i < m && !NUM(f); i++ ) {    for ( i = 1, f = tps[0]; i < m && !NUM(f); i++ ) {
                 uezgcdpz(vl,f,tps[i],&t); f = t;      uezgcdpz(vl,f,tps[i],&t); f = t;
         }    }
         *r = f;    *r = f;
 }  }
   
 void gcdmonomial(vl,dcl,m,r)  void gcdmonomial(vl,dcl,m,r)
Line 140  DCP *dcl;
Line 191  DCP *dcl;
 int m;  int m;
 P *r;  P *r;
 {  {
         int i,j,n;    int i,j,n;
         P g,x,s,t;    P g,x,s,t;
         Q d;    Q d;
         DCP dc;    DCP dc;
         VN vn;    VN vn;
   
         for ( i = 0; i < m; i++ )    for ( i = 0; i < m; i++ )
                 if ( !dcl[i] ) {      if ( !dcl[i] ) {
                         *r = (P)ONE; return;        *r = (P)ONE; return;
                 }      }
         for ( n = 0, dc = dcl[0]; dc; dc = NEXT(dc), n++ );    for ( n = 0, dc = dcl[0]; dc; dc = NEXT(dc), n++ );
         vn = (VN)ALLOCA(n*sizeof(struct oVN));    vn = (VN)ALLOCA(n*sizeof(struct oVN));
         for ( i = 0, dc = dcl[0]; i < n; dc = NEXT(dc), i++ ) {    for ( i = 0, dc = dcl[0]; i < n; dc = NEXT(dc), i++ ) {
                 vn[i].v = VR(COEF(dc)); vn[i].n = QTOS(DEG(dc));      vn[i].v = VR(COEF(dc)); vn[i].n = QTOS(DEG(dc));
         }    }
         for ( i = 1; i < m; i++ ) {    for ( i = 1; i < m; i++ ) {
                 for ( j = 0; j < n; j++ ) {      for ( j = 0; j < n; j++ ) {
                         for ( dc = dcl[i]; dc; dc = NEXT(dc) )        for ( dc = dcl[i]; dc; dc = NEXT(dc) )
                                 if ( VR(COEF(dc)) == vn[j].v ) {          if ( VR(COEF(dc)) == vn[j].v ) {
                                         vn[j].n = MIN(vn[j].n,QTOS(DEG(dc))); break;            vn[j].n = MIN(vn[j].n,QTOS(DEG(dc))); break;
                                 }          }
                         if ( !dc )        if ( !dc )
                                 vn[j].n = 0;          vn[j].n = 0;
                 }      }
                 for ( j = n-1; j >= 0 && !vn[j].n; j-- );      for ( j = n-1; j >= 0 && !vn[j].n; j-- );
                 if ( j < 0 ) {      if ( j < 0 ) {
                         *r = (P)ONE; return;        *r = (P)ONE; return;
                 } else      } else
                         n = j+1;        n = j+1;
         }    }
         for ( j = 0, g = (P)ONE; j < n; j++ )    for ( j = 0, g = (P)ONE; j < n; j++ )
                 if ( vn[j].n ) {      if ( vn[j].n ) {
                         MKV(vn[j].v,x); STOQ(vn[j].n,d);        MKV(vn[j].v,x); STOQ(vn[j].n,d);
                         pwrp(vl,x,d,&t); mulp(vl,t,g,&s); g = s;        pwrp(vl,x,d,&t); mulp(vl,t,g,&s); g = s;
                 }      }
         *r = g;    *r = g;
 }  }
   
 void nezgcdnpzmain(vl,pl,m,pr)  void nezgcdnpzmain(vl,pl,m,pr)
Line 183  VL vl;
Line 234  VL vl;
 P *pl,*pr;  P *pl,*pr;
 int m;  int m;
 {  {
         P *ppl,*pcl;    P *ppl,*pcl;
         int i;    int i;
         P cont,gcd,t,s;    P cont,gcd,t,s;
         DCP dc;    DCP dc;
   
         ppl = (P *)ALLOCA(m*sizeof(P));    ppl = (P *)ALLOCA(m*sizeof(P));
         pcl = (P *)ALLOCA(m*sizeof(P));    pcl = (P *)ALLOCA(m*sizeof(P));
         for ( i = 0; i < m; i++ )    for ( i = 0; i < m; i++ )
                 pcp(vl,pl[i],&ppl[i],&pcl[i]);      pcp(vl,pl[i],&ppl[i],&pcl[i]);
         nezgcdnpz(vl,pcl,m,&cont);    nezgcdnpz(vl,pcl,m,&cont);
         sqfrp(vl,ppl[0],&dc);    sqfrp(vl,ppl[0],&dc);
         for ( dc = NEXT(dc), gcd = (P)ONE; dc; dc = NEXT(dc) ) {    for ( dc = NEXT(dc), gcd = (P)ONE; dc; dc = NEXT(dc) ) {
                 if ( NUM(COEF(dc)) )      if ( NUM(COEF(dc)) )
                         continue;        continue;
                 nezgcdnpp(vl,dc,ppl+1,m-1,&t);      nezgcdnpp(vl,dc,ppl+1,m-1,&t);
                 if ( NUM(t) )      if ( NUM(t) )
                         continue;        continue;
                 mulp(vl,gcd,t,&s); gcd = s;      mulp(vl,gcd,t,&s); gcd = s;
         }    }
         mulp(vl,gcd,cont,pr);    mulp(vl,gcd,cont,pr);
 }  }
   
 void nezgcdnpp(vl,dc,pl,m,r)  void nezgcdnpp(vl,dc,pl,m,r)
Line 212  P *pl;
Line 263  P *pl;
 int m;  int m;
 P *r;  P *r;
 {  {
         int i,k;    int i,k;
         P g,t,s,gcd;    P g,t,s,gcd;
         P *pm;    P *pm;
   
         nezgcdnp_sqfr_primitive(vl,COEF(dc),pl,m,&gcd);    nezgcdnp_sqfr_primitive(vl,COEF(dc),pl,m,&gcd);
         if ( NUM(gcd) ) {    if ( NUM(gcd) ) {
                 *r = (P)ONE; return;      *r = (P)ONE; return;
         }    }
         pm = (P *) ALLOCA(m*sizeof(P));    pm = (P *) ALLOCA(m*sizeof(P));
         for ( i = 0; i < m; i++ ) {    for ( i = 0; i < m; i++ ) {
                 divsp(vl,pl[i],gcd,&pm[i]);      divsp(vl,pl[i],gcd,&pm[i]);
                 if ( NUM(pm[i]) ) {      if ( NUM(pm[i]) ) {
                         *r = gcd; return;        *r = gcd; return;
                 }      }
         }    }
         for ( g = gcd, k = QTOS(DEG(dc))-1; k > 0; k-- ) {    for ( g = gcd, k = QTOS(DEG(dc))-1; k > 0; k-- ) {
                 nezgcdnp_sqfr_primitive(vl,g,pm,m,&t);      nezgcdnp_sqfr_primitive(vl,g,pm,m,&t);
                 if ( NUM(t) )      if ( NUM(t) )
                         break;        break;
                 mulp(vl,gcd,t,&s); gcd = s;      mulp(vl,gcd,t,&s); gcd = s;
                 for ( i = 0; i < m; i++ ) {      for ( i = 0; i < m; i++ ) {
                         divsp(vl,pm[i],t,&s);        divsp(vl,pm[i],t,&s);
                         if ( NUM(s) ) {        if ( NUM(s) ) {
                                 *r = gcd; return;          *r = gcd; return;
                         }        }
                         pm[i] = s;        pm[i] = s;
                 }      }
         }    }
         *r = gcd;    *r = gcd;
 }  }
   
 /*  /*
Line 256  VL vl;
Line 307  VL vl;
 int m;  int m;
 P p0,*ps,*pr;  P p0,*ps,*pr;
 {  {
         /* variables */    /* variables */
         P p00,g,h,g0,h0,a0,b0;    P p00,g,h,g0,h0,a0,b0;
         P lp0,lgp0,lp00,lg,lg0,cbd,tq,t;    P lp0,lgp0,lp00,lg,lg0,cbd,tq,t;
         P *lc;    P *lc;
         P *tps;    P *tps;
         VL nvl1,nvl2,nvl,tvl;    VL nvl1,nvl2,nvl,tvl;
         V v;    V v;
         int i,j,k,d0,dg,dg0,dmin;    int i,j,k,d0,dg,dg0,dmin,z;
         VN vn0,vn1,vnt;    VN vn1;
         int nv,flag;    int nv,flag,max;
   
         /* set-up */    /* set-up */
         if ( NUM(p0) ) {    if ( NUM(p0) ) {
                 *pr = (P) ONE; return;      *pr = (P) ONE; return;
         }    }
         for ( v = VR(p0), i = 0; i < m; i++ )    for ( v = VR(p0), i = 0; i < m; i++ )
                 if ( NUM(ps[i]) || (v != VR(ps[i])) ) {      if ( NUM(ps[i]) || (v != VR(ps[i])) ) {
                         *pr = (P)ONE; return;        *pr = (P)ONE; return;
                 }      }
         tps = (P *) ALLOCA(m*sizeof(P));    tps = (P *) ALLOCA(m*sizeof(P));
         for ( i = 0; i < m; i++ )    for ( i = 0; i < m; i++ )
                 tps[i] = ps[i];      tps[i] = ps[i];
         sortplist(tps,m);    sortplist(tps,m);
         /* deg(tps[0]) <= deg(tps[1]) <= ... */    /* deg(tps[0]) <= deg(tps[1]) <= ... */
   
         if ( !cmpq(DEG(DC(p0)),ONE) ) {    if ( !cmpq(DEG(DC(p0)),ONE) ) {
                 if ( divcheck(vl,tps,m,(P)ONE,p0) )      if ( divcheck(vl,tps,m,(P)ONE,p0) )
                         *pr = p0;        *pr = p0;
                 else      else
                         *pr = (P)ONE;        *pr = (P)ONE;
                 return;      return;
         }    }
   
         lp0 = LC(p0); dmin = d0 = deg(v,p0); lc = (P *)ALLOCA((m+1)*sizeof(P));    lp0 = LC(p0); dmin = d0 = deg(v,p0); lc = (P *)ALLOCA((m+1)*sizeof(P));
         for ( lc[0] = lp0, i = 0; i < m; i++ )    for ( lc[0] = lp0, i = 0; i < m; i++ )
                 lc[i+1] = LC(tps[i]);      lc[i+1] = LC(tps[i]);
         clctv(vl,p0,&nvl);    clctv(vl,p0,&nvl);
         for ( i = 0; i < m; i++ ) {    for ( i = 0; i < m; i++ ) {
                 clctv(vl,tps[i],&nvl1); mergev(vl,nvl,nvl1,&nvl2); nvl = nvl2;      clctv(vl,tps[i],&nvl1); mergev(vl,nvl,nvl1,&nvl2); nvl = nvl2;
         }    }
         nezgcdnpz(nvl,lc,m+1,&lg);    nezgcdnpz(nvl,lc,m+1,&lg);
   
         mulp(nvl,p0,lg,&lgp0); k = dbound(v,lgp0)+1; cbound(nvl,lgp0,(Q *)&cbd);    mulp(nvl,p0,lg,&lgp0); k = dbound(v,lgp0)+1; cbound(nvl,lgp0,(Q *)&cbd);
         for ( nv = 0, tvl = nvl; tvl; tvl = NEXT(tvl), nv++ );    for ( nv = 0, tvl = nvl; tvl; tvl = NEXT(tvl), nv++ );
         W_CALLOC(nv,struct oVN,vn0); W_CALLOC(nv,struct oVN,vnt);    W_CALLOC(nv,struct oVN,vn1);
         W_CALLOC(nv,struct oVN,vn1);    for ( i = 0, tvl = NEXT(nvl); tvl; tvl = NEXT(tvl), i++ )
         for ( i = 0, tvl = NEXT(nvl); tvl; tvl = NEXT(tvl), i++ )      vn1[i].v = tvl->v;
                 vn1[i].v = vn0[i].v = tvl->v;  
   
         /* main loop */    /* main loop */
         for ( dg = deg(v,tps[0]) + 1; ; next(vn0) )    /* finally, 'max' random evaluations will be generated */
                 do {    for ( max = 1, dg = deg(v,tps[0]) + 1; ; max = 2*max )
                         for ( i = 0, j = 0; vn0[i].v; i++ )      for ( z = 0; z < max; z++ ) {
                                 if ( vn0[i].n ) vnt[j++].v = (V)i;        for ( i = 0; vn1[i].v; i++ )
                         vnt[j].n = 0;          vn1[i].n = mt_genrand()%max;
   
                         /* find b s.t. LC(p0)(b), LC(tps[i])(b) != 0 */        /* find b s.t. LC(p0)(b), LC(tps[i])(b) != 0 */
                         mulsgn(vn0,vnt,j,vn1); substvp(nvl,p0,vn1,&p00);  
                         flag = (!zerovpchk(nvl,lp0,vn1) && sqfrchk(p00));        substvp(nvl,p0,vn1,&p00);
                         for ( i = 0; flag && (i < m); i++ )        flag = (!zerovpchk(nvl,lp0,vn1) && sqfrchk(p00));
                                 flag &= (!zerovpchk(nvl,LC(tps[i]),vn1));        for ( i = 0; flag && (i < m); i++ )
                         if ( !flag )          flag &= (!zerovpchk(nvl,LC(tps[i]),vn1));
                                 continue;        if ( !flag )
           continue;
   
                         /* substitute y -> b */        /* substitute y -> b */
                         substvp(nvl,lg,vn1,&lg0); lp00 = LC(p00);        substvp(nvl,lg,vn1,&lg0); lp00 = LC(p00);
                         /* extended-gcd in 1-variable */        /* extended-gcd in 1-variable */
                         uexgcdnp(nvl,p00,tps,m,vn1,(Q)cbd,&g0,&h0,&a0,&b0,(Q *)&tq);        uexgcdnp(nvl,p00,tps,m,vn1,(Q)cbd,&g0,&h0,&a0,&b0,(Q *)&tq);
                         if ( NUM(g0) ) {        if ( NUM(g0) ) {
                                 *pr = (P)ONE; return;          *pr = (P)ONE; return;
                         } else if ( dg > ( dg0 = deg(v,g0) ) ) {        } else if ( dg > ( dg0 = deg(v,g0) ) ) {
                                 dg = dg0;          dg = dg0;
                                 if ( dg == d0 ) {          if ( dg == d0 ) {
                                         if ( divcheck(nvl,tps,m,lp0,p0) ) {            if ( divcheck(nvl,tps,m,lp0,p0) ) {
                                                 *pr = p0; return;              *pr = p0; return;
                                         }            }
                                 } else if ( dg == deg(v,tps[0]) ) {          } else if ( dg == deg(v,tps[0]) ) {
                                         if ( divtpz(nvl,p0,tps[0],&t) &&            if ( divtpz(nvl,p0,tps[0],&t) &&
                                                 divcheck(nvl,tps+1,m-1,LC(tps[0]),tps[0]) ) {              divcheck(nvl,tps+1,m-1,LC(tps[0]),tps[0]) ) {
                                                 *pr = tps[0]; return;              *pr = tps[0]; return;
                                         } else            } else
                                                 break;              break;
                                 } else {          } else {
                                         henmv(nvl,vn1,lgp0,g0,h0,a0,b0,lg,lp0,lg0,lp00,(Q)tq,k,&g,&h);            henmv(nvl,vn1,lgp0,g0,h0,a0,b0,lg,lp0,lg0,lp00,(Q)tq,k,&g,&h);
                                         if ( divtpz(nvl,lgp0,g,&t) &&            if ( divtpz(nvl,lgp0,g,&t) &&
                                                 divcheck(nvl,tps,m,lg,g) ) {              divcheck(nvl,tps,m,lg,g) ) {
                                                 pcp(nvl,g,pr,&t); return;              pcp(nvl,g,pr,&t); return;
                                         }            }
                                 }          }
                         }        }
                 } while ( !nextbin(vnt,j) );      }
 }  }
   
 void intersectv(vl1,vl2,vlp)  void intersectv(vl1,vl2,vlp)
 VL vl1,vl2,*vlp;  VL vl1,vl2,*vlp;
 {  {
         int i,n;    int i,n;
         VL tvl;    VL tvl;
         VN tvn;    VN tvn;
   
         if ( !vl1 || !vl2 ) {    if ( !vl1 || !vl2 ) {
                 *vlp = 0; return;      *vlp = 0; return;
         }    }
         for ( n = 0, tvl = vl1; tvl; tvl = NEXT(tvl), n++ );    for ( n = 0, tvl = vl1; tvl; tvl = NEXT(tvl), n++ );
         tvn = (VN) ALLOCA(n*sizeof(struct oVN));    tvn = (VN) ALLOCA(n*sizeof(struct oVN));
         for ( i = 0, tvl = vl1; i < n; tvl = NEXT(tvl), i++ ) {    for ( i = 0, tvl = vl1; i < n; tvl = NEXT(tvl), i++ ) {
                 tvn[i].v = tvl->v; tvn[i].n = 0;      tvn[i].v = tvl->v; tvn[i].n = 0;
         }    }
         for ( tvl = vl2; tvl; tvl = NEXT(tvl) )    for ( tvl = vl2; tvl; tvl = NEXT(tvl) )
                 for ( i = 0; i < n; i++ )      for ( i = 0; i < n; i++ )
                         if ( tvn[i].v == tvl->v ) {        if ( tvn[i].v == tvl->v ) {
                                 tvn[i].n = 1; break;          tvn[i].n = 1; break;
                         }        }
         vntovl(tvn,n,vlp);    vntovl(tvn,n,vlp);
 }  }
   
 int pcoef(vl,ivl,p,cp)  int pcoef(vl,ivl,p,cp)
Line 377  VL vl,ivl;
Line 428  VL vl,ivl;
 P p;  P p;
 P *cp;  P *cp;
 {  {
         VL nvl,tvl,svl,mvl,mvl0;    VL nvl,tvl,svl,mvl,mvl0;
         P t;    P t;
   
         if ( NUM(p) ) {    if ( NUM(p) ) {
                 cp[0] = p; return 1;      cp[0] = p; return 1;
         } else {    } else {
                 clctv(vl,p,&nvl);      clctv(vl,p,&nvl);
                 for ( tvl = nvl, mvl0 = 0; tvl; tvl = NEXT(tvl) ) {      for ( tvl = nvl, mvl0 = 0; tvl; tvl = NEXT(tvl) ) {
                         for ( svl = ivl; svl; svl = NEXT(svl) )        for ( svl = ivl; svl; svl = NEXT(svl) )
                                 if ( tvl->v == svl->v )          if ( tvl->v == svl->v )
                                         break;            break;
                         if ( !svl ) {        if ( !svl ) {
                                 if ( !mvl0 ) {          if ( !mvl0 ) {
                                         NEWVL(mvl0); mvl = mvl0;            NEWVL(mvl0); mvl = mvl0;
                                 } else {          } else {
                                         NEWVL(NEXT(mvl)); mvl = NEXT(mvl);            NEWVL(NEXT(mvl)); mvl = NEXT(mvl);
                                 }          }
                                 mvl->v = tvl->v;          mvl->v = tvl->v;
                         }        }
                 }      }
                 if ( mvl0 )      if ( mvl0 )
                         NEXT(mvl) = ivl;        NEXT(mvl) = ivl;
                 else      else
                         mvl0 = ivl;        mvl0 = ivl;
                 reorderp(mvl0,nvl,p,&t);      reorderp(mvl0,nvl,p,&t);
                 return pcoef0(mvl0,ivl,t,cp);      return pcoef0(mvl0,ivl,t,cp);
         }    }
 }  }
   
 int pcoef0(vl,ivl,p,cp)  int pcoef0(vl,ivl,p,cp)
Line 411  VL vl,ivl;
Line 462  VL vl,ivl;
 P p;  P p;
 P *cp;  P *cp;
 {  {
         int cn,n;    int cn,n;
         DCP dc;    DCP dc;
         V v;    V v;
         VL tvl;    VL tvl;
   
         if ( NUM(p) ) {    if ( NUM(p) ) {
                 cp[0] = p; return 1;      cp[0] = p; return 1;
         } else {    } else {
                 for ( v = VR(p), tvl = ivl; tvl; tvl = NEXT(tvl) )      for ( v = VR(p), tvl = ivl; tvl; tvl = NEXT(tvl) )
                         if ( v == tvl->v )        if ( v == tvl->v )
                                 break;          break;
                 if ( tvl ) {      if ( tvl ) {
                         cp[0] = p; return 1;        cp[0] = p; return 1;
                 } else {      } else {
                         for ( dc = DC(p), n = 0; dc; dc = NEXT(dc) ) {        for ( dc = DC(p), n = 0; dc; dc = NEXT(dc) ) {
                                 cn = pcoef0(vl,ivl,COEF(dc),cp); cp += cn; n += cn;          cn = pcoef0(vl,ivl,COEF(dc),cp); cp += cn; n += cn;
                         }        }
                         return n;        return n;
                 }      }
         }    }
 }  }
   
 int lengthp(p)  int lengthp(p)
 P p;  P p;
 {  {
         int n;    int n;
         DCP dc;    DCP dc;
   
         if ( NUM(p) )    if ( NUM(p) )
                 return 1;      return 1;
         else {    else {
                 for ( dc = DC(p), n = 0; dc; dc = NEXT(dc) )      for ( dc = DC(p), n = 0; dc; dc = NEXT(dc) )
                         n += lengthp(COEF(dc));        n += lengthp(COEF(dc));
                 return n;      return n;
         }    }
   }
   
   int nonzero_const_term(P p)
   {
     DCP dc;
   
     if ( !p )
       return 0;
     else if ( NUM(p) )
       return 1;
     else {
       dc = DC(p);
       for ( ; NEXT(dc); dc = NEXT(dc) );
       if ( DEG(dc) )
         return 0;
       else
         return nonzero_const_term(COEF(dc));
     }
 }  }

Legend:
Removed from v.1.1  
changed lines
  Added in v.1.7

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