version 1.10, 2002/11/22 07:32:10 |
version 1.16, 2003/01/06 09:23:27 |
|
|
/* $OpenXM: OpenXM_contrib2/asir2000/engine/Fgfs.c,v 1.9 2002/11/01 06:47:41 noro Exp $ */ |
/* $OpenXM: OpenXM_contrib2/asir2000/engine/Fgfs.c,v 1.15 2003/01/04 09:06:17 noro Exp $ */ |
|
|
#include "ca.h" |
#include "ca.h" |
|
|
Line 12 void partial_sqfrsf(VL vl,V v,P f,P *r,DCP *dcp); |
|
Line 12 void partial_sqfrsf(VL vl,V v,P f,P *r,DCP *dcp); |
|
void gcdsf(VL vl,P *pa,int k,P *r); |
void gcdsf(VL vl,P *pa,int k,P *r); |
void mfctrsfmain(VL vl, P f, DCP *dcp); |
void mfctrsfmain(VL vl, P f, DCP *dcp); |
void next_evaluation_point(int *mev,int n); |
void next_evaluation_point(int *mev,int n); |
void estimatelc_sf(VL vl,VL rvl,P c,DCP dc,P *lcp); |
void estimatelc_sf(VL vl,VL rvl,P c,DCP dc,int *mev,P *lcp); |
void mfctrsf_hensel(VL vl,VL rvl,P f,P pp0,P u0,P v0,P lcu,P lcv,P *up); |
void mfctrsf_hensel(VL vl,VL rvl,P f,P pp0,P u0,P v0,P lcu,P lcv,int *mev,P *up); |
void substvp_sf(VL vl,VL rvl,P f,int *mev,P *r); |
void substvp_sf(VL vl,VL rvl,P f,int *mev,P *r); |
void shift_sf(VL vl, VL rvl, P f, int *mev, int sgn, P *r); |
void shift_sf(VL vl, VL rvl, P f, int *mev, int sgn, P *r); |
void adjust_coef_sf(VL vl,VL rvl,P lcu,P u0,P *r); |
void adjust_coef_sf(VL vl,VL rvl,P lcu,P u0,int *mev,P *r); |
void extended_gcd_modyk(P u0,P v0,V x,V y,int dy,P *cu,P *cv); |
void extended_gcd_modyk(P u0,P v0,V x,V y,int dy,P *cu,P *cv); |
void poly_to_gfsn_poly(VL vl,P f,V v,P *r); |
void poly_to_gfsn_poly(VL vl,P f,V v,P *r); |
void gfsn_poly_to_poly(VL vl,P f,V v,P *r); |
void gfsn_poly_to_poly(VL vl,P f,V v,P *r); |
Line 50 void sqfrsf(VL vl, P f, DCP *dcp) |
|
Line 50 void sqfrsf(VL vl, P f, DCP *dcp) |
|
{ |
{ |
DCP dc,dct; |
DCP dc,dct; |
Obj obj; |
Obj obj; |
P t,s,c; |
P t,s,c,cont; |
VL tvl,nvl; |
VL tvl,onevl; |
|
|
simp_ff((Obj)f,&obj); f = (P)obj; |
simp_ff((Obj)f,&obj); f = (P)obj; |
lex_lc(f,&c); divsp(vl,f,c,&t); f = t; |
lex_lc(f,&c); divsp(vl,f,c,&t); f = t; |
monomialfctr(vl,f,&t,&dc); f = t; |
monomialfctr(vl,f,&t,&dc); f = t; |
clctv(vl,f,&tvl); vl = tvl; |
clctv(vl,f,&tvl); vl = tvl; |
|
NEWVL(onevl); NEXT(onevl)=0; |
if ( !vl ) |
if ( !vl ) |
; |
; |
else if ( !NEXT(vl) ) { |
else if ( !NEXT(vl) ) { |
Line 65 void sqfrsf(VL vl, P f, DCP *dcp) |
|
Line 66 void sqfrsf(VL vl, P f, DCP *dcp) |
|
} else { |
} else { |
t = f; |
t = f; |
for ( tvl = vl; tvl; tvl = NEXT(tvl) ) { |
for ( tvl = vl; tvl; tvl = NEXT(tvl) ) { |
reordvar(vl,tvl->v,&nvl); |
onevl->v = tvl->v; |
cont_pp_mv_sf(vl,NEXT(nvl),t,&c,&s); t = s; |
cont_pp_mv_sf(vl,onevl,t,&cont,&s); t = s; |
sqfrsf(vl,c,&dct); |
sqfrsf(vl,cont,&dct); |
dc = append_dc(dc,NEXT(dct)); |
dc = append_dc(dc,NEXT(dct)); |
} |
} |
sqfrsfmain(vl,t,&dct); |
sqfrsfmain(vl,t,&dct); |
Line 473 void cont_pp_mv_sf(VL vl,VL rvl,P p,P *c,P *pp) |
|
Line 474 void cont_pp_mv_sf(VL vl,VL rvl,P p,P *c,P *pp) |
|
MP t; |
MP t; |
int i,m; |
int i,m; |
P *ps; |
P *ps; |
|
struct order_spec spec, currentspec; |
|
extern struct order_spec dp_current_spec; |
|
|
|
currentspec = dp_current_spec; |
|
create_order_spec(0,&spec); |
|
initd(&spec); |
ptod(vl,rvl,p,&dp); |
ptod(vl,rvl,p,&dp); |
for ( t = BDY(dp), m = 0; t; t = NEXT(t), m++ ); |
for ( t = BDY(dp), m = 0; t; t = NEXT(t), m++ ); |
ps = (P *)ALLOCA(m*sizeof(P)); |
ps = (P *)ALLOCA(m*sizeof(P)); |
Line 481 void cont_pp_mv_sf(VL vl,VL rvl,P p,P *c,P *pp) |
|
Line 487 void cont_pp_mv_sf(VL vl,VL rvl,P p,P *c,P *pp) |
|
ps[i] = C(t); |
ps[i] = C(t); |
gcdsf(vl,ps,m,c); |
gcdsf(vl,ps,m,c); |
divsp(vl,p,*c,pp); |
divsp(vl,p,*c,pp); |
|
initd(¤tspec); |
} |
} |
|
|
void mfctrsf(VL vl, P f, DCP *dcp) |
void mfctrsf(VL vl, P f, DCP *dcp) |
Line 515 void mfctrsfmain(VL vl, P f, DCP *dcp) |
|
Line 522 void mfctrsfmain(VL vl, P f, DCP *dcp) |
|
P gcd,g,df,dfmin; |
P gcd,g,df,dfmin; |
P pa[2]; |
P pa[2]; |
P f0,pp0,spp0,c,c0,x,y,u,v,lcf,lcu,lcv,u0,v0,t,s; |
P f0,pp0,spp0,c,c0,x,y,u,v,lcf,lcu,lcv,u0,v0,t,s; |
P ype,yme; |
P ype,yme,fin; |
GFS ev,evy; |
GFS ev,evy; |
P *fp0; |
P *fp0; |
int *mev,*win; |
int *mev,*win; |
Line 556 void mfctrsfmain(VL vl, P f, DCP *dcp) |
|
Line 563 void mfctrsfmain(VL vl, P f, DCP *dcp) |
|
} |
} |
pa[0] = f; |
pa[0] = f; |
pa[1] = dfmin; |
pa[1] = dfmin; |
gcdsf_main(vl,pa,2,&gcd); |
gcdsf(vl,pa,2,&gcd); |
if ( !NUM(gcd) ) { |
if ( !NUM(gcd) ) { |
/* f = gcd * f/gcd */ |
/* f = gcd * f/gcd */ |
mfctrsfmain(vl,gcd,&dc1); |
mfctrsfmain(vl,gcd,&dc1); |
Line 578 void mfctrsfmain(VL vl, P f, DCP *dcp) |
|
Line 585 void mfctrsfmain(VL vl, P f, DCP *dcp) |
|
} |
} |
NEXT(tvl) = 0; |
NEXT(tvl) = 0; |
|
|
|
fin = f; |
reorderp(nvl,vl,f,&g); f = g; |
reorderp(nvl,vl,f,&g); f = g; |
vx = nvl->v; |
vx = nvl->v; |
vy = NEXT(nvl)->v; |
vy = NEXT(nvl)->v; |
Line 626 void mfctrsfmain(VL vl, P f, DCP *dcp) |
|
Line 634 void mfctrsfmain(VL vl, P f, DCP *dcp) |
|
|
|
if ( !NEXT(dc) ) { |
if ( !NEXT(dc) ) { |
/* f is irreducible */ |
/* f is irreducible */ |
NEWDC(dc); DEG(dc) = ONE; COEF(dc) = f; NEXT(dc) = 0; |
NEWDC(dc); DEG(dc) = ONE; COEF(dc) = fin; NEXT(dc) = 0; |
*dcp = dc; |
*dcp = dc; |
return; |
return; |
} |
} |
Line 645 void mfctrsfmain(VL vl, P f, DCP *dcp) |
|
Line 653 void mfctrsfmain(VL vl, P f, DCP *dcp) |
|
lcf = COEF(DC(f)); |
lcf = COEF(DC(f)); |
mfctrsf(nvl,lcf,&dct); |
mfctrsf(nvl,lcf,&dct); |
/* skip the first element (= a number) */ |
/* skip the first element (= a number) */ |
dct = NEXT(dct); |
lcfdc = NEXT(dct); |
|
|
/* shift lcfdc; c <- c(X+mev) */ |
|
for ( lcfdc = 0; dct; dct = NEXT(dct) ) { |
|
NEXTDC(lcfdc,dcs); |
|
DEG(dcs) = DEG(dct); |
|
shift_sf(nvl,rvl,COEF(dct),mev,1,&COEF(dcs)); |
|
} |
|
NEXT(dcs) = 0; |
|
|
|
/* np = number of bivariate factors */ |
/* np = number of bivariate factors */ |
for ( np = 0, dct = dc; dct; dct = NEXT(dct), np++ ); |
for ( np = 0, dct = dc; dct; dct = NEXT(dct), np++ ); |
fp0 = (P *)ALLOCA((np+1)*sizeof(P)); |
fp0 = (P *)ALLOCA((np+1)*sizeof(P)); |
Line 664 void mfctrsfmain(VL vl, P f, DCP *dcp) |
|
Line 664 void mfctrsfmain(VL vl, P f, DCP *dcp) |
|
l = tl = (P *)ALLOCA((np+1)*sizeof(P)); |
l = tl = (P *)ALLOCA((np+1)*sizeof(P)); |
win = W_ALLOC(np+1); |
win = W_ALLOC(np+1); |
|
|
/* f <- f(X+mev) */ |
|
shift_sf(nvl,rvl,f,mev,1,&s); f = s; |
|
|
|
for ( k = 1, win[0] = 1, --np; ; ) { |
for ( k = 1, win[0] = 1, --np; ; ) { |
itogfs(1,&u0); |
itogfs(1,&u0); |
/* u0 = product of selected factors */ |
/* u0 = product of selected factors */ |
Line 675 void mfctrsfmain(VL vl, P f, DCP *dcp) |
|
Line 672 void mfctrsfmain(VL vl, P f, DCP *dcp) |
|
} |
} |
/* we have to consider the content */ |
/* we have to consider the content */ |
/* f0 = c0*u0*v0 */ |
/* f0 = c0*u0*v0 */ |
mulp(nvl,LC(u0),c0,&c); estimatelc_sf(nvl,rvl,c,lcfdc,&lcu); |
mulp(nvl,LC(u0),c0,&c); estimatelc_sf(nvl,rvl,c,lcfdc,mev,&lcu); |
divsp(nvl,pp0,u0,&v0); |
divsp(nvl,pp0,u0,&v0); |
mulp(nvl,LC(v0),c0,&c); estimatelc_sf(nvl,rvl,c,lcfdc,&lcv); |
mulp(nvl,LC(v0),c0,&c); estimatelc_sf(nvl,rvl,c,lcfdc,mev,&lcv); |
mfctrsf_hensel(nvl,rvl,f,pp0,u0,v0,lcu,lcv,&u); |
mfctrsf_hensel(nvl,rvl,f,pp0,u0,v0,lcu,lcv,mev,&u); |
if ( u ) { |
if ( u ) { |
/* save the factor */ |
/* save the factor */ |
reorderp(vl,nvl,u,&t); |
reorderp(vl,nvl,u,&t); |
/* x -> x-mev, y -> y-evy */ |
/* y -> y-evy */ |
shift_sf(vl,rvl,t,mev,-1,&s); substp(vl,s,vy,yme,tl++); |
substp(vl,t,vy,yme,tl++); |
|
|
/* update f,pp0 */ |
/* update f,pp0 */ |
divsp(nvl,f,u,&t); f = t; |
divsp(nvl,f,u,&t); f = t; |
Line 714 void mfctrsfmain(VL vl, P f, DCP *dcp) |
|
Line 711 void mfctrsfmain(VL vl, P f, DCP *dcp) |
|
} |
} |
} |
} |
reorderp(vl,nvl,f,&t); |
reorderp(vl,nvl,f,&t); |
/* x -> x-mev, y -> y-evy */ |
/* y -> y-evy */ |
shift_sf(vl,rvl,t,mev,-1,&s); substp(vl,s,vy,yme,tl++); |
substp(vl,t,vy,yme,tl++); |
*tl = 0; |
*tl = 0; |
for ( dc0 = 0, i = 0; l[i]; i++ ) { |
for ( dc0 = 0, i = 0; l[i]; i++ ) { |
NEXTDC(dc0,dc); DEG(dc) = ONE; COEF(dc) = l[i]; |
NEXTDC(dc0,dc); DEG(dc) = ONE; COEF(dc) = l[i]; |
Line 742 void next_evaluation_point(int *e,int n) |
|
Line 739 void next_evaluation_point(int *e,int n) |
|
|
|
/* |
/* |
* dc : f1^E1*...*fk^Ek |
* dc : f1^E1*...*fk^Ek |
* find e1,...,ek s.t. fi(0)^ei | c |
* find e1,...,ek s.t. fi(mev)^ei | c |
* and return f1^e1*...*fk^ek |
* and return f1^e1*...*fk^ek |
* vl = (vx,vy,rvl) |
* vl = (vx,vy,rvl) |
*/ |
*/ |
|
|
void estimatelc_sf(VL vl,VL rvl,P c,DCP dc,P *lcp) |
void estimatelc_sf(VL vl,VL rvl,P c,DCP dc,int *mev,P *lcp) |
{ |
{ |
DCP dct; |
DCP dct; |
P r,c1,c2,t,s,f; |
P r,c1,c2,t,s,f; |
Line 758 void estimatelc_sf(VL vl,VL rvl,P c,DCP dc,P *lcp) |
|
Line 755 void estimatelc_sf(VL vl,VL rvl,P c,DCP dc,P *lcp) |
|
if ( NUM(COEF(dct)) ) |
if ( NUM(COEF(dct)) ) |
continue; |
continue; |
/* constant part */ |
/* constant part */ |
substvp_sf(vl,rvl,COEF(dct),0,&f); |
substvp_sf(vl,rvl,COEF(dct),mev,&f); |
d = QTOS(DEG(dct)); |
d = QTOS(DEG(dct)); |
for ( i = 0, c1 = c; i < d; i++ ) |
for ( i = 0, c1 = c; i < d; i++ ) |
if ( !divtp(vl,c1,f,&c2) ) |
if ( !divtp(vl,c1,f,&c2) ) |
Line 819 void shift_sf(VL vl, VL rvl, P f, int *mev, int sgn, P |
|
Line 816 void shift_sf(VL vl, VL rvl, P f, int *mev, int sgn, P |
|
* pp(f(0)) = u0*v0 |
* pp(f(0)) = u0*v0 |
*/ |
*/ |
|
|
void mfctrsf_hensel(VL vl,VL rvl,P f,P pp0,P u0,P v0,P lcu,P lcv,P *up) |
void mfctrsf_hensel(VL vl,VL rvl,P f,P pp0,P u0,P v0,P lcu,P lcv,int *mev,P *up) |
{ |
{ |
VL tvl,onevl; |
VL tvl,onevl; |
P t,s,w,u,v,ff,si,wu,wv,fj,cont; |
P t,s,w,u,v,ff,si,wu,wv,fj,cont; |
Line 828 void mfctrsf_hensel(VL vl,VL rvl,P f,P pp0,P u0,P v0,P |
|
Line 825 void mfctrsf_hensel(VL vl,VL rvl,P f,P pp0,P u0,P v0,P |
|
int dy,n,i,dbd,nv,j; |
int dy,n,i,dbd,nv,j; |
int *md; |
int *md; |
P *uh,*vh; |
P *uh,*vh; |
P x,du0,dv0,m,q,r; |
P x,du0,dv0,m,q,r,fin; |
P *cu,*cv; |
P *cu,*cv; |
GFSN inv; |
GFSN inv; |
|
|
/* adjust coeffs */ |
/* check the validity of lc's and adjust coeffs */ |
/* u0 = am x^m+ ... -> lcu*x^m + a(m-1)*(lcu(0)/am)*x^(m-1)+... */ |
|
/* v0 = bm x^l+ ... -> lcv*x^l + b(l-1)*(lcv(0)/bl)*x^(l-1)+... */ |
|
/* f -> lcu*lcv*x^(m+l)+... */ |
/* f -> lcu*lcv*x^(m+l)+... */ |
adjust_coef_sf(vl,rvl,lcu,u0,&u); |
mulp(vl,lcu,lcv,&t); |
adjust_coef_sf(vl,rvl,lcv,v0,&v); |
if ( !divtp(vl,t,LC(f),&m) ) { |
mulp(vl,lcu,lcv,&t); divsp(vl,t,LC(f),&m); mulp(vl,m,f,&t); f = t; |
*up = 0; return; |
|
} |
|
mulp(vl,m,f,&t); f = t; |
|
/* u0 = am x^m+ ... -> lcu*x^m + a(m-1)*(lcu(mev)/am)*x^(m-1)+... */ |
|
/* v0 = bm x^l+ ... -> lcv*x^l + b(l-1)*(lcv(mev)/bl)*x^(l-1)+... */ |
|
adjust_coef_sf(vl,rvl,lcu,u0,mev,&u); |
|
adjust_coef_sf(vl,rvl,lcv,v0,mev,&v); |
|
|
|
/* f <- f(X+mev), u <- u(X+mev), v <- v(X+mev) */ |
|
fin = f; |
|
shift_sf(vl,rvl,f,mev,1,&s); f = s; |
|
shift_sf(vl,rvl,u,mev,1,&s); u = s; |
|
shift_sf(vl,rvl,v,mev,1,&s); v = s; |
|
|
vx = vl->v; vy = NEXT(vl)->v; |
vx = vl->v; vy = NEXT(vl)->v; |
n = getdeg(vx,f); |
n = getdeg(vx,f); |
dy = getdeg(vy,f)+1; |
dy = getdeg(vy,f)+1; |
Line 930 void mfctrsf_hensel(VL vl,VL rvl,P f,P pp0,P u0,P v0,P |
|
Line 937 void mfctrsf_hensel(VL vl,VL rvl,P f,P pp0,P u0,P v0,P |
|
} |
} |
if ( !wu ) { |
if ( !wu ) { |
gfsn_poly_to_poly(vl,u,vy,&t); |
gfsn_poly_to_poly(vl,u,vy,&t); |
if ( divtp(vl,f,t,&q) ) { |
shift_sf(vl,rvl,t,mev,-1,&s); |
cont_pp_mv_sf(vl,onevl,t,&cont,up); |
if ( divtp(vl,fin,s,&q) ) { |
|
cont_pp_mv_sf(vl,onevl,s,&cont,up); |
return; |
return; |
} |
} |
} |
} |
if ( !wv ) { |
if ( !wv ) { |
gfsn_poly_to_poly(vl,v,vy,&t); |
gfsn_poly_to_poly(vl,v,vy,&t); |
if ( divtp(vl,f,t,&q) ) { |
shift_sf(vl,rvl,t,mev,-1,&s); |
|
if ( divtp(vl,fin,s,&q) ) { |
cont_pp_mv_sf(vl,onevl,q,&cont,up); |
cont_pp_mv_sf(vl,onevl,q,&cont,up); |
return; |
return; |
} |
} |
Line 949 void mfctrsf_hensel(VL vl,VL rvl,P f,P pp0,P u0,P v0,P |
|
Line 958 void mfctrsf_hensel(VL vl,VL rvl,P f,P pp0,P u0,P v0,P |
|
} |
} |
} |
} |
|
|
void adjust_coef_sf(VL vl,VL rvl,P lcu,P u0,P *r) |
void adjust_coef_sf(VL vl,VL rvl,P lcu,P u0,int *mev,P *r) |
{ |
{ |
P lcu0,cu; |
P lcu0,cu; |
DCP dc0,dcu,dc; |
DCP dc0,dcu,dc; |
|
|
substvp_sf(vl,rvl,lcu,0,&lcu0); |
substvp_sf(vl,rvl,lcu,mev,&lcu0); |
divsp(vl,lcu0,LC(u0),&cu); |
divsp(vl,lcu0,LC(u0),&cu); |
for ( dc0 = 0, dcu = DC(u0); dcu; dcu = NEXT(dcu) ) { |
for ( dc0 = 0, dcu = DC(u0); dcu; dcu = NEXT(dcu) ) { |
if ( !dc0 ) { |
if ( !dc0 ) { |