=================================================================== RCS file: /home/cvs/OpenXM_contrib2/asir2000/engine/Fgfs.c,v retrieving revision 1.2 retrieving revision 1.3 diff -u -p -r1.2 -r1.3 --- OpenXM_contrib2/asir2000/engine/Fgfs.c 2002/09/26 09:07:42 1.2 +++ OpenXM_contrib2/asir2000/engine/Fgfs.c 2002/09/27 08:40:48 1.3 @@ -1,11 +1,16 @@ -/* $OpenXM: OpenXM_contrib2/asir2000/engine/Fgfs.c,v 1.1 2002/09/26 04:33:16 noro Exp $ */ +/* $OpenXM: OpenXM_contrib2/asir2000/engine/Fgfs.c,v 1.2 2002/09/26 09:07:42 noro Exp $ */ #include "ca.h" +void cont_pp_mv_sf(VL vl,VL rvl,P p,P *c,P *pp); +void gcdsf_main(VL vl,P *pa,int m,P *r); +void ugcdsf(P *pa,int m,P *r); +void head_monomial(V v,P p,P *coef,P *term); void gcdsf(VL vl,P *pa,int k,P *r) { - P *ps,*pl,*pm,*cp; + P *ps,*pl,*pm; + P **cp; int *cn; DCP *ml; Obj obj; @@ -17,21 +22,21 @@ void gcdsf(VL vl,P *pa,int k,P *r) for ( i = 0, m = 0; i < k; i++ ) { simp_ff((Obj)pa[i],&obj); if ( obj ) - ps[m++] = obj; + ps[m++] = (P)obj; } if ( !m ) { *r = 0; return; } if ( m == 1 ) { - *r = BDY(n0); + *r = ps[0]; return; } pl = (P *)ALLOCA(m*sizeof(P)); ml = (DCP *)ALLOCA(m*sizeof(DCP)); for ( i = 0; i < m; i++ ) monomialfctr(vl,ps[i],&pl[i],&ml[i]); - gcdmonomial(vl,ml,m,&mg); simp_ff(mg,&mgsf); + gcdmonomial(vl,ml,m,&mg); simp_ff((Obj)mg,&obj); mgsf = (P)obj; for ( i = 0, nvl = vl, avl = 0; nvl && i < m; i++ ) { clctv(vl,pl[i],&tvl); intersectv(nvl,tvl,&svl); nvl = svl; @@ -73,7 +78,7 @@ void gcdsf(VL vl,P *pa,int k,P *r) void ugcdsf(P *pa,int m,P *r) { - P *pa; + P *ps; int i; UM w1,w2,w3,w; int d; @@ -83,10 +88,15 @@ void ugcdsf(P *pa,int m,P *r) *r = pa[0]; return; } + for ( i = 0; i < m; i++ ) + if ( NUM(pa[i]) ) { + itogfs(1,r); + return; + } ps = (P *)ALLOCA(m*sizeof(P)); - v = VR(pa[0]); sort_by_deg(m,pa,ps); - d = getdeg(ps[m-1],v); + v = VR(ps[m-1]); + d = getdeg(v,ps[m-1]); w1 = W_UMALLOC(d); w2 = W_UMALLOC(d); w3 = W_UMALLOC(d); @@ -96,7 +106,7 @@ void ugcdsf(P *pa,int m,P *r) gcdsfum(w1,w2,w3); w = w1; w1 = w3; w3 = w; if ( !DEG(w1) ) { - MKGFS(0,*r); + itogfs(1,r); return; } } @@ -108,10 +118,13 @@ void ugcdsf(P *pa,int m,P *r) void gcdsf_main(VL vl,P *pa,int m,P *r) { - int nv,i,i0,imin,d,d0,d1,d2; - V v0,vmin; + int nv,i,i0,imin,d,d0,d1,d2,dmin,index; + V v,v0,vmin; VL tvl,nvl,rvl,nvl0,rvl0; - P *pc, *ps, *ph; + P *pc, *ps, *ph,*lps; + P x,t,cont,hg,g,hm,mod,s; + P hge,ge,ce,he,u,cof1e,mode,mod1,adj,cof1,coadj,q; + GFS sf; for ( nv = 0, tvl = vl; tvl; tvl = NEXT(tvl), nv++); if ( nv == 1 ) { @@ -124,7 +137,7 @@ void gcdsf_main(VL vl,P *pa,int m,P *r) v = tvl->v; i = 0; do { - d = getdeg(pa[i],v); + d = getdeg(v,pa[i]); if ( i == 0 || (d < d0) ) { d0 = d; i0 = i; v0 = v; } @@ -141,9 +154,9 @@ void gcdsf_main(VL vl,P *pa,int m,P *r) NEXTVL(rvl0,rvl); rvl->v = tvl->v; } /* rvl = remaining variables */ - NEXT(rvl) = 0; + NEXT(rvl) = 0; rvl = rvl0; /* nvl = ...,vmin */ - NEXTVL(nvl0,nvl); nvl->v = vmin; NEXT(nvl) = 0; + NEXTVL(nvl0,nvl); nvl->v = vmin; NEXT(nvl) = 0; nvl = nvl0; MKV(vmin,x); /* for content and primitive part */ @@ -152,7 +165,7 @@ void gcdsf_main(VL vl,P *pa,int m,P *r) ph = (P *)ALLOCA(m*sizeof(P)); /* separate the contents */ for ( i = 0; i < m; i++ ) { - reorderp(vl,nvl,pa[i],&t); + reorderp(nvl,vl,pa[i],&t); cont_pp_mv_sf(nvl,rvl,t,&pc[i],&ps[i]); head_monomial(vmin,ps[i],&ph[i],&t); } @@ -166,13 +179,13 @@ void gcdsf_main(VL vl,P *pa,int m,P *r) while ( 1 ) { g = 0; - cofmin = 0; + cof1 = 0; hm = 0; - MKGFS(0,mod); + itogfs(1,&mod); index = 0; - while ( getdeg(mod,vmin) <= d+1 ) { + for ( index = 0; getdeg(vmin,mod) <= d+1; index++ ) { /* evaluation pt */ - MKGFS(index,s); + indextogfs(index,&s); substp(nvl,hg,vmin,s,&hge); if ( !hge ) continue; @@ -181,42 +194,43 @@ void gcdsf_main(VL vl,P *pa,int m,P *r) /* ge = GCD(ps[0]|x=s,...,ps[m-1]|x=s) */ gcdsf(nvl,ph,m,&ge); head_monomial(vmin,ge,&ce,&he); - if ( NUM(hge) ) { + if ( NUM(he) ) { *r = cont; return; } - divgfs(hge,ce,&t); mulp(nvl,t,ge,&u); ge = u; + divgfs((GFS)hge,(GFS)ce,&sf); t = (P)sf; + mulp(nvl,t,ge,&u); ge = u; divsp(nvl,ph[imin],ge,&t); mulp(nvl,hge,t,&cof1e); /* hm=0 : reset; he==hm : lucky */ - if ( !hm || !cmpp(he,hm) ) { + if ( !hm || !compp(nvl,he,hm) ) { substp(nvl,mod,vmin,s,&mode); divsp(nvl,mod,mode,&mod1); /* adj = mod/(mod|x=s)*(ge-g|x=s) */ substp(nvl,g,vmin,s,&t); subp(nvl,ge,t,&u); mulp(nvl,mod1,u,&adj); /* coadj = mod/(mod|vmin=s)*(cof1e-cof1e|vmin=s) */ substp(nvl,cof1,vmin,s,&t); - subp(nvl,cof1,t,&u); mulp(nvl,mod1,u,&coadj); + subp(nvl,cof1e,t,&u); mulp(nvl,mod1,u,&coadj); if ( !adj ) { /* adj == gcd ? */ for ( i = 0; i < m; i++ ) - if ( !divtp(nvl,lps[i],adj,&t) ) + if ( !divtp(nvl,lps[i],g,&t) ) break; if ( i == m ) { - cont_pp_mv_sf(nvl,rvl,adj,&t,&u); + cont_pp_mv_sf(nvl,rvl,g,&t,&u); mulp(nvl,cont,u,&t); - reorderp(nvl,vl,t,r); + reorderp(vl,nvl,t,r); return; } } else if ( !coadj ) { - /* ps[vmin]/coadj == gcd ? */ - if ( divtp(nvl,lps[vmin],coadj,&q) ) { + /* ps[imin]/coadj == gcd ? */ + if ( divtp(nvl,lps[imin],cof1,&q) ) { for ( i = 0; i < m; i++ ) if ( !divtp(nvl,lps[i],q,&t) ) break; if ( i == m ) { cont_pp_mv_sf(nvl,rvl,q,&t,&u); mulp(nvl,cont,u,&t); - reorderp(nvl,vl,t,r); + reorderp(vl,nvl,t,r); return; } } @@ -239,7 +253,7 @@ void gcdsf_main(VL vl,P *pa,int m,P *r) /* d1==d2, but hm!=he => both are unlucky */ g = 0; cof1 = 0; - MKGFS(0,mod); + itogfs(1,&mod); } } } @@ -251,16 +265,18 @@ void head_monomial(V v,P p,P *coef,P *term) P t,s,u; DCP dc; GFS one; + VL vl; - MKGFS(0,one); t = one; + itogfs(1,&one); + t = (P)one; while ( 1 ) { if ( NUM(p) || VR(p) == v ) { *coef = p; *term = t; return; } else { - NEWDC(dc); MKGFS(0,one); - COEF(dc) = one; DEG(dc) = DEG(DC(p)); + NEWDC(dc); + COEF(dc) = (P)one; DEG(dc) = DEG(DC(p)); MKP(VR(p),dc,s); mulp(vl,t,s,&u); t = u; p = COEF(DC(p)); @@ -278,8 +294,8 @@ void cont_pp_mv_sf(VL vl,VL rvl,P p,P *c,P *pp) ptod(vl,rvl,p,&dp); for ( t = BDY(dp), m = 0; t; t = NEXT(t), m++ ); ps = (P *)ALLOCA(m*sizeof(P)); - for ( t = BDY(dp), i = 0; t; t = NEXT(t), i++ ); + for ( t = BDY(dp), i = 0; t; t = NEXT(t), i++ ) ps[i] = C(t); ugcdsf(ps,m,c); - divsp(vl,p,c,pp); + divsp(vl,p,*c,pp); }