version 1.17, 2001/09/03 07:01:06 |
version 1.27, 2002/11/01 05:43:35 |
|
|
/* $OpenXM: OpenXM_contrib2/asir2000/engine/Hgfs.c,v 1.16 2001/09/03 01:04:26 noro Exp $ */ |
/* $OpenXM: OpenXM_contrib2/asir2000/engine/Hgfs.c,v 1.26 2002/10/25 02:43:40 noro Exp $ */ |
|
|
#include "ca.h" |
#include "ca.h" |
|
#include "inline.h" |
|
|
struct p_pair { |
void lnfsf(int n,UM p0,UM p1,struct p_pair *list,UM np0,UM np1); |
UM p0; |
void extractcoefbm(BM f,int dx,UM r); |
UM p1; |
|
struct p_pair *next; |
|
}; |
|
|
|
void canzassf(UM,int,UM *); |
int comp_dum(DUM a,DUM b) |
void lnfsf(int,UM,UM,struct p_pair *,UM,UM); |
|
void minipolysf(UM,UM,UM); |
|
void czsfum(UM,UM *); |
|
void gensqfrsfum(UM,DUM); |
|
void sfbmtop(BM,V,V,P *); |
|
void cont_pp_sfp(VL,P,P *,P *); |
|
void sfcsump(VL,P,P *); |
|
void mulsfbmarray(int,BM,ML,int,int *,V,V,P *); |
|
void const_term(P,UM); |
|
|
|
void sfbsqfr(P,V,V,DCP *); |
|
void sfusqfr(P,DCP *); |
|
|
|
int comp_dum(a,b) |
|
DUM a,b; |
|
{ |
{ |
if ( DEG(a->f) > DEG(b->f) ) |
if ( DEG(a->f) > DEG(b->f) ) |
return -1; |
return -1; |
|
|
return 0; |
return 0; |
} |
} |
|
|
void fctrsf(p,dcp) |
void ufctrsf(P p,DCP *dcp) |
P p; |
|
DCP *dcp; |
|
{ |
{ |
int n,i,j,k; |
int n,i,j,k; |
DCP dc,dc0; |
DCP dc,dc0; |
P lc; |
P lc; |
P zp; |
|
UM mp; |
UM mp; |
UM *tl; |
UM *tl; |
|
Obj obj; |
struct oDUM *udc,*udc1; |
struct oDUM *udc,*udc1; |
|
|
simp_ff(p,&zp); p = zp; |
simp_ff((Obj)p,&obj); p = (P)obj; |
if ( !p ) { |
if ( !p ) { |
*dcp = 0; return; |
NEWDC(dc); COEF(dc) = 0; DEG(dc) = ONE; |
|
NEXT(dc) = 0; *dcp = dc; |
|
return; |
} |
} |
mp = W_UMALLOC(UDEG(p)); |
mp = W_UMALLOC(UDEG(p)); |
ptosfum(p,mp); |
ptosfum(p,mp); |
|
|
NEXT(dc) = 0; *dcp = dc0; |
NEXT(dc) = 0; *dcp = dc0; |
} |
} |
|
|
void gensqfrsfum(p,dc) |
void gensqfrsfum(UM p,struct oDUM *dc) |
UM p; |
|
struct oDUM *dc; |
|
{ |
{ |
int n,i,j,d,mod; |
int n,i,j,d,mod; |
UM t,s,g,f,f1,b; |
UM t,s,g,f,f1,b; |
Line 147 struct oDUM *dc; |
|
Line 128 struct oDUM *dc; |
|
} |
} |
} |
} |
|
|
void randsfum(d,p) |
void randsfum(int d,UM p) |
int d; |
|
UM p; |
|
{ |
{ |
int i; |
int i; |
|
|
|
|
p->d = i; |
p->d = i; |
} |
} |
|
|
void pwrmodsfum(p,e,f,pr) |
void pwrmodsfum(UM p,int e,UM f,UM pr) |
int e; |
|
UM p,f,pr; |
|
{ |
{ |
UM wt,ws,q; |
UM wt,ws,q; |
|
|
|
|
} |
} |
} |
} |
|
|
void spwrsfum(m,f,e,r) |
void spwrsfum(UM m,UM f,N e,UM r) |
UM f,m,r; |
|
N e; |
|
{ |
{ |
UM t,s,q; |
UM t,s,q; |
N e1; |
N e1; |
|
|
} |
} |
} |
} |
|
|
void tracemodsfum(m,f,e,r) |
void tracemodsfum(UM m,UM f,int e,UM r) |
UM f,m,r; |
|
int e; |
|
{ |
{ |
UM t,s,q,u; |
UM t,s,q,u; |
int i; |
int i; |
|
|
cpyum(s,r); |
cpyum(s,r); |
} |
} |
|
|
void make_qmatsf(p,tab,mp) |
void make_qmatsf(UM p,UM *tab,int ***mp) |
UM p; |
|
UM *tab; |
|
int ***mp; |
|
{ |
{ |
int n,i,j; |
int n,i,j; |
int *c; |
int *c; |
|
|
mat[i][i] = _subsf(mat[i][i],one); |
mat[i][i] = _subsf(mat[i][i],one); |
} |
} |
|
|
void nullsf(mat,n,ind) |
void nullsf(int **mat,int n,int *ind) |
int **mat; |
|
int *ind; |
|
int n; |
|
{ |
{ |
int i,j,l,s,h,inv; |
int i,j,l,s,h,inv; |
int *t,*u; |
int *t,*u; |
|
|
} |
} |
} |
} |
|
|
void null_to_solsf(mat,ind,n,r) |
void null_to_solsf(int **mat,int *ind,int n,UM *r) |
int **mat; |
|
int *ind; |
|
int n; |
|
UM *r; |
|
{ |
{ |
int i,j,k,l; |
int i,j,k,l; |
int *c; |
int *c; |
Line 327 nullsf(mat,n,ind) |
|
Line 290 nullsf(mat,n,ind) |
|
null_to_solsf(ind,n,r) |
null_to_solsf(ind,n,r) |
*/ |
*/ |
|
|
void czsfum(f,r) |
void czsfum(UM f,UM *r) |
UM f,*r; |
|
{ |
{ |
int i,j; |
int i,j; |
int d,n,ord; |
int d,n,ord; |
|
|
r[j] = 0; |
r[j] = 0; |
} |
} |
|
|
int berlekampsf(p,df,tab,r) |
int berlekampsf(UM p,int df,UM *tab,UM *r) |
UM p; |
|
int df; |
|
UM *tab,*r; |
|
{ |
{ |
int n,i,j,k,nf,d,nr; |
int n,i,j,k,nf,d,nr; |
int **mat; |
int **mat; |
|
|
} |
} |
} |
} |
} |
} |
|
/* NOT REACHED */ |
|
error("berlekampsf : cannot happen"); |
|
return 0; |
} |
} |
|
|
void minipolysf(f,p,mp) |
void minipolysf(UM f,UM p,UM mp) |
UM f,p,mp; |
|
{ |
{ |
struct p_pair *list,*l,*l1,*lprev; |
struct p_pair *list,*l,*l1,*lprev; |
int n,d; |
int n,d; |
|
|
} |
} |
} |
} |
|
|
void lnfsf(n,p0,p1,list,np0,np1) |
void lnfsf(int n,UM p0,UM p1,struct p_pair *list,UM np0,UM np1) |
int n; |
|
UM p0,p1; |
|
struct p_pair *list; |
|
UM np0,np1; |
|
{ |
{ |
int inv,h,d1; |
int h,d1; |
UM t0,t1,s0,s1; |
UM t0,t1,s0,s1; |
struct p_pair *l; |
struct p_pair *l; |
|
|
|
|
} |
} |
} |
} |
|
|
int find_rootsf(p,root) |
int find_rootsf(UM p,int *root) |
UM p; |
|
int *root; |
|
{ |
{ |
UM *r; |
UM *r; |
int i,j,n; |
int i,n; |
|
|
n = DEG(p); |
n = DEG(p); |
r = ALLOCA((DEG(p))*sizeof(UM)); |
r = ALLOCA((DEG(p))*sizeof(UM)); |
|
|
return n; |
return n; |
} |
} |
|
|
void canzassf(f,d,r) |
void canzassf(UM f,int d,UM *r) |
UM f,*r; |
|
int d; |
|
{ |
{ |
UM t,s,u,w,g,o; |
UM t,s,u,w,g,o; |
N n1,n2,n3,n4,n5; |
N n1,n2,n3,n4,n5; |
UM *b; |
UM *b; |
int n,m,i,q,ed; |
int n,q,ed; |
|
|
if ( DEG(f) == d ) { |
if ( DEG(f) == d ) { |
r[0] = UMALLOC(d); cpyum(f,r[0]); |
r[0] = UMALLOC(d); cpyum(f,r[0]); |
Line 560 void ptosfbm(int,P,BM); |
|
Line 513 void ptosfbm(int,P,BM); |
|
|
|
/* f = f(x,y) */ |
/* f = f(x,y) */ |
|
|
void sfhensel(count,f,x,evp,sfp,listp) |
void sfhensel(int count,P f,V x,int degbound,GFS *evp,P *sfp,ML *listp) |
int count; |
|
P f; |
|
V x; |
|
GFS *evp; |
|
P *sfp; |
|
ML *listp; |
|
{ |
{ |
int i,j; |
int i; |
int fn,n; |
int fn; |
ML rlist; |
ML rlist; |
BM fl; |
BM fl; |
VL vl,nvl; |
VL vl,nvl; |
V y; |
V y; |
int dx,dy,mev; |
int dx,dy,bound; |
GFS ev; |
GFS ev; |
P f1,t,yev,c,sf; |
P f1,t,c,sf; |
DCP dc; |
DCP dc,dct,dc0; |
UM w,w1,q,fm,hm; |
UM q,fm,hm; |
UM *gm; |
UM *gm; |
struct oEGT tmp0,tmp1,eg_hensel,eg_hensel_t; |
struct oEGT tmp0,tmp1,eg_hensel,eg_hensel_t; |
|
|
|
|
*listp = rlist = MLALLOC(1); rlist->n = fn; rlist->c[0] = 0; |
*listp = rlist = MLALLOC(1); rlist->n = fn; rlist->c[0] = 0; |
return; |
return; |
} |
} |
/* pass the the leading coeff. to the first element */ |
if ( degbound >= 0 ) { |
c = dc->c; dc = NEXT(dc); |
/* |
mulp(vl,dc->c,c,&t); dc->c = t; |
* reconstruct dc so that |
|
* dc[1],... : factors satisfy degree bound |
|
* dc[0] : product of others |
|
*/ |
|
c = dc->c; dc = NEXT(dc); |
|
dc0 = 0; |
|
fn = 0; |
|
while ( dc ) { |
|
if ( getdeg(x,COEF(dc)) <= degbound ) { |
|
dct = NEXT(dc); NEXT(dc) = dc0; dc0 = dc; dc = dct; |
|
fn++; |
|
} else { |
|
mulp(vl,COEF(dc),c,&t); c = t; |
|
dc = NEXT(dc); |
|
} |
|
} |
|
if ( OID(c) == O_P ) { |
|
NEWDC(dc); COEF(dc) = c; DEG(dc) = ONE; NEXT(dc) = dc0; |
|
fn++; |
|
} else { |
|
mulp(vl,dc0->c,c,&t); dc0->c = t; dc = dc0; |
|
} |
|
} else { |
|
/* pass the the leading coeff. to the first element */ |
|
c = dc->c; dc = NEXT(dc); |
|
mulp(vl,dc->c,c,&t); dc->c = t; |
|
} |
|
|
/* convert mod y-a factors into UM */ |
/* convert mod y-a factors into UM */ |
gm = (UM *)ALLOCA(fn*sizeof(UM)); |
gm = (UM *)ALLOCA(fn*sizeof(UM)); |
|
|
ptosfum(dc->c,gm[i]); |
ptosfum(dc->c,gm[i]); |
} |
} |
|
|
|
/* set bound */ |
|
/* g | f, lc_y(g) = lc_y(f) => deg_y(g) <= deg_y(f) */ |
|
/* so, bound = dy is sufficient, but we use slightly large value */ |
|
bound = dy+2; |
|
|
/* f(x,y) -> f(x,y+ev) */ |
/* f(x,y) -> f(x,y+ev) */ |
fl = BMALLOC(dx,dy); |
fl = BMALLOC(dx,bound); |
ptosfbm(dy,f,fl); |
ptosfbm(bound,f,fl); |
if ( ev ) shiftsfbm(fl,FTOIF(CONT(ev))); |
if ( ev ) shiftsfbm(fl,FTOIF(CONT(ev))); |
|
|
/* sf = f(x+ev) */ |
/* sf = f(x+ev) */ |
|
|
hm = W_UMALLOC(dx); |
hm = W_UMALLOC(dx); |
|
|
q = W_UMALLOC(dx); |
q = W_UMALLOC(dx); |
rlist = MLALLOC(fn); rlist->n = fn; rlist->bound = dy; |
rlist = MLALLOC(fn); rlist->n = fn; rlist->bound = bound; |
fprintf(asir_out,"%d candidates\n",fn); |
fprintf(asir_out,"%d candidates\n",fn); |
init_eg(&eg_hensel); |
init_eg(&eg_hensel); |
for ( i = 0; i < fn-1; i++ ) { |
for ( i = 0; i < fn-1; i++ ) { |
|
|
get_eg(&tmp0); |
get_eg(&tmp0); |
/* fl = gm[i]*hm mod y */ |
/* fl = gm[i]*hm mod y */ |
divsfum(fm,gm[i],hm); |
divsfum(fm,gm[i],hm); |
/* fl is replaced by the cofactor of gk mod y^dy */ |
/* fl is replaced by the cofactor of gk mod y^bound */ |
/* rlist->c[i] = gk */ |
/* rlist->c[i] = gk */ |
sfhenmain2(fl,gm[i],hm,dy,(BM *)&rlist->c[i]); |
sfhenmain2(fl,gm[i],hm,bound,(BM *)&rlist->c[i]); |
cpyum(hm,fm); |
cpyum(hm,fm); |
get_eg(&tmp1); add_eg(&eg_hensel_t,&tmp0,&tmp1); |
get_eg(&tmp1); add_eg(&eg_hensel_t,&tmp0,&tmp1); |
add_eg(&eg_hensel,&tmp0,&tmp1); |
add_eg(&eg_hensel,&tmp0,&tmp1); |
|
|
|
|
/* main variable of f = x */ |
/* main variable of f = x */ |
|
|
int sfberle(vl,f,count,ev,dcp) |
int sfberle(VL vl,P f,int count,GFS *ev,DCP *dcp) |
VL vl; |
|
P f; |
|
int count; |
|
GFS *ev; |
|
DCP *dcp; |
|
{ |
{ |
UM wf,wf1,wf2,wfs,gcd; |
UM wf,wf1,wf2,wfs,gcd; |
ML flist; |
int fn,n; |
int fn,fn1,n; |
|
GFS m,fm; |
GFS m,fm; |
DCP dc,dct,dc0; |
DCP dc,dct,dc0; |
VL nvl; |
VL nvl; |
V x,y; |
V x,y; |
P g,lc,lc0,f0; |
P lc,lc0,f0; |
int j,q1,index,i; |
Obj obj; |
|
int j,q,index,i; |
|
|
clctv(vl,f,&nvl); vl = nvl; |
clctv(vl,f,&nvl); vl = nvl; |
x = vl->v; y = vl->next->v; |
x = vl->v; y = vl->next->v; |
simp_ff(f,&g); g = f; |
simp_ff((Obj)f,&obj); f = (P)obj; |
n = QTOS(DEG(DC(f))); |
n = QTOS(DEG(DC(f))); |
wf = W_UMALLOC(n); wf1 = W_UMALLOC(n); wf2 = W_UMALLOC(n); |
wf = W_UMALLOC(n); wf1 = W_UMALLOC(n); wf2 = W_UMALLOC(n); |
wfs = W_UMALLOC(n); gcd = W_UMALLOC(n); |
wfs = W_UMALLOC(n); gcd = W_UMALLOC(n); |
q1 = field_order_sf()-1; |
q = field_order_sf(); |
lc = DC(f)->c; |
lc = DC(f)->c; |
for ( j = 0, fn = n + 1, index = 0; |
for ( j = 0, fn = n + 1, index = 0; |
index < q1 && j < count && fn > 1; index++ ) { |
index < q && j < count && fn > 1; index++ ) { |
MKGFS(index,m); |
indextogfs(index,&m); |
substp(vl,lc,y,(P)m,&lc0); |
substp(vl,lc,y,(P)m,&lc0); |
if ( lc0 ) { |
if ( lc0 ) { |
substp(vl,f,y,(P)m,&f0); |
substp(vl,f,y,(P)m,&f0); |
ptosfum(f0,wf); cpyum(wf,wf1); |
ptosfum(f0,wf); cpyum(wf,wf1); |
diffsfum(wf1,wf2); gcdsfum(wf1,wf2,gcd); |
diffsfum(wf1,wf2); gcdsfum(wf1,wf2,gcd); |
if ( DEG(gcd) == 0 ) { |
if ( DEG(gcd) == 0 ) { |
fctrsf(f0,&dc); |
ufctrsf(f0,&dc); |
for ( dct = NEXT(dc), i = 0; dct; dct = NEXT(dct), i++ ); |
for ( dct = NEXT(dc), i = 0; dct; dct = NEXT(dct), i++ ); |
if ( i < fn ) { |
if ( i < fn ) { |
dc0 = dc; fn = i; fm = m; |
dc0 = dc; fn = i; fm = m; |
|
|
} |
} |
} |
} |
} |
} |
if ( index == q1 ) |
if ( index == q ) |
return 0; |
return 0; |
else if ( fn == 1 ) |
else if ( fn == 1 ) |
return 1; |
return 1; |
|
|
} |
} |
} |
} |
|
|
void sfgcdgen(f,blist,clistp) |
void sfgcdgen(P f,ML blist,ML *clistp) |
P f; |
|
ML blist,*clistp; |
|
{ |
{ |
int i; |
int i; |
int n,d,np; |
int n,d,np; |
Line 747 ML blist,*clistp; |
|
Line 718 ML blist,*clistp; |
|
|
|
/* f = g0*h0 mod y -> f = gk*hk mod y^(dy+1), f is replaced by hk */ |
/* f = g0*h0 mod y -> f = gk*hk mod y^(dy+1), f is replaced by hk */ |
|
|
void sfhenmain2(f,g0,h0,dy,gp) |
void sfhenmain2(BM f,UM g0,UM h0,int dy,BM *gp) |
BM f; |
|
UM g0,h0; |
|
int dy; |
|
BM *gp; |
|
{ |
{ |
int i,j,k,l; |
int i,k; |
int *px,*py; |
int dx; |
int **pp,**pp1; |
UM wt,wa,wb,q,w1,w2,wh1,wg1,ws; |
int dx,np,dr,tmp; |
|
UM wt,wa,wb,wm,q,w1,w2,wh1,wg1,ws; |
|
UM wc,wd,we,wz; |
UM wc,wd,we,wz; |
BM wb0,wb1; |
BM wb0,wb1; |
int dg,dh; |
int dg,dh; |
|
|
cpyum(COEF(hk)[i],COEF(f)[i]); |
cpyum(COEF(hk)[i],COEF(f)[i]); |
} |
} |
|
|
|
/* a0*g+b0*h = 1 mod y -> a*g+b*h = 1 mod y^(dy+1) */ |
|
|
|
void sfexgcd_by_hensel(BM g,BM h,int dy,BM *ap,BM *bp) |
|
{ |
|
int i,k; |
|
int dx; |
|
UM wt,wa,wb,q,w1,w2,ws; |
|
UM wc,wd,we,wz,wa1,wb1; |
|
BM wz0,wz1; |
|
int dg,dh; |
|
BM a,b,c; |
|
|
|
dg = degbm(g); |
|
dh = degbm(h); |
|
dx = dg+dh; |
|
|
|
a = BMALLOC(dh,dy); |
|
b = BMALLOC(dg,dy); |
|
/* c holds a*g+b*h-1 */ |
|
c = BMALLOC(dg+dh,dy); |
|
|
|
W_BMALLOC(dx,dy,wz0); W_BMALLOC(dx,dy,wz1); |
|
|
|
wt = W_UMALLOC(dx); ws = W_UMALLOC(dx); q = W_UMALLOC(2*dx); |
|
wa1 = W_UMALLOC(2*dx); wb1 = W_UMALLOC(2*dx); |
|
wc = W_UMALLOC(2*dx); wd = W_UMALLOC(2*dx); |
|
we = W_UMALLOC(2*dx); wz = W_UMALLOC(2*dx); |
|
|
|
/* compute wa,wb s.t. wa*g0+wb*h0 = 1 mod y */ |
|
w1 = W_UMALLOC(dg); cpyum(COEF(g)[0],w1); |
|
w2 = W_UMALLOC(dh); cpyum(COEF(h)[0],w2); |
|
wa = W_UMALLOC(2*dx); wb = W_UMALLOC(2*dx); /* XXX */ |
|
eucsfum(w1,w2,wa,wb); |
|
cpyum(wa,COEF(a)[0]); cpyum(wb,COEF(b)[0]); |
|
|
|
/* initialize c to a*g+b*h-1 */ |
|
mulsfbm(a,g,c); mulsfbm(b,h,wz0); addtosfbm(wz0,c); |
|
COEF(COEF(c)[0])[0] = 0; |
|
|
|
fprintf(stderr,"dy=%d\n",dy); |
|
for ( k = 1; k <= dy; k++ ) { |
|
fprintf(stderr,"."); |
|
|
|
/* at this point, a*g+b*h = 1 mod y^k, c = a*g+b*h-1 */ |
|
|
|
/* wt = -((a*g+b*h-1)/y^k) */ |
|
cpyum(COEF(c)[k],wt); |
|
for ( i = DEG(wt); i >= 0; i-- ) |
|
COEF(wt)[i] = _chsgnsf(COEF(wt)[i]); |
|
|
|
/* compute wa1,wb1 s.t. wa1*g0+wb1*h0 = wt */ |
|
mulsfum(wa,wt,wa1); DEG(wa1) = divsfum(wa1,COEF(h)[0],q); |
|
mulsfum(wa1,COEF(g)[0],wc); subsfum(wt,wc,wd); |
|
DEG(wd) = divsfum(wd,COEF(h)[0],wb1); |
|
|
|
/* c += ((wa1*g+wb1*h)*y^k mod y^(dy+1) */ |
|
/* wz0 = wa1*y^k */ |
|
clearbm(dx,wz0); |
|
cpyum(wa1,COEF(wz0)[k]); |
|
|
|
/* wz1 = wz0*g mod y^(dy+1) */ |
|
clearbm(dx,wz1); |
|
mulsfbm(g,wz0,wz1); |
|
/* c += wz1 */ |
|
addtosfbm(wz1,c); |
|
|
|
/* wz0 = wb1*y^k */ |
|
clearbm(dx,wz0); |
|
cpyum(wb1,COEF(wz0)[k]); |
|
|
|
/* wz1 = wz0*h mod y^(dy+1) */ |
|
clearbm(dx,wz1); |
|
mulsfbm(h,wz0,wz1); |
|
/* c += wz1 */ |
|
addtosfbm(wz1,c); |
|
|
|
/* a += wa1*y^k, b += wb1*y^k */ |
|
cpyum(wa1,COEF(a)[k]); |
|
cpyum(wb1,COEF(b)[k]); |
|
} |
|
fprintf(stderr,"\n"); |
|
DEG(a) = dy; |
|
DEG(b) = dy; |
|
*ap = a; |
|
*bp = b; |
|
} |
|
|
/* fl->c[i] = coef_y(f,i) */ |
/* fl->c[i] = coef_y(f,i) */ |
|
|
void ptosfbm(dy,f,fl) |
void ptosfbm(int dy,P f,BM fl) |
int dy; |
|
P f; |
|
BM fl; |
|
{ |
{ |
DCP dc; |
DCP dc; |
int d,i,dx; |
int d,i,dx; |
|
|
|
|
/* x : main variable */ |
/* x : main variable */ |
|
|
void sfbmtop(f,x,y,fp) |
void sfbmtop(BM f,V x,V y,P *fp) |
BM f; |
|
V x,y; |
|
P *fp; |
|
{ |
{ |
UM *c; |
UM *c; |
int i,j,d,a,dy; |
int i,j,d,a,dy; |
|
|
if ( DEG(c[j]) >= i && (a = COEF(c[j])[i]) ) { |
if ( DEG(c[j]) >= i && (a = COEF(c[j])[i]) ) { |
NEWDC(dct); |
NEWDC(dct); |
STOQ(j,DEG(dct)); |
STOQ(j,DEG(dct)); |
MKGFS(IFTOF(a),b); |
iftogfs(a,&b); |
COEF(dct) = (P)b; |
COEF(dct) = (P)b; |
NEXT(dct) = dc; |
NEXT(dct) = dc; |
dc = dct; |
dc = dct; |
|
|
*fp = 0; |
*fp = 0; |
} |
} |
|
|
void sfsqfr(f,dcp) |
void sfsqfr(P f,DCP *dcp) |
P f; |
|
DCP *dcp; |
|
{ |
{ |
P g; |
Obj obj; |
V x; |
|
DCP dc; |
DCP dc; |
VL vl; |
VL vl; |
|
|
simp_ff(f,&g); f = g; |
simp_ff((Obj)f,&obj); f = (P)obj; |
clctv(CO,f,&vl); |
clctv(CO,f,&vl); |
if ( !vl ) { |
if ( !vl ) { |
/* f is a const */ |
/* f is a const */ |
NEWDC(dc); DEG(dc) = ONE; COEF(dc) = f; NEXT(dc) = 0; *dcp = dc; |
NEWDC(dc); DEG(dc) = ONE; COEF(dc) = f; NEXT(dc) = 0; *dcp = dc; |
} else if ( !NEXT(vl) ) |
} else if ( !NEXT(vl) ) |
sfusqfr(f,dcp); |
sfusqfr(f,dcp); |
else if ( !NEXT(NEXT(vl)) ) |
|
sfbsqfr(f,x,NEXT(vl)->v,dcp); |
|
else |
else |
error("sfsqfr : not implemented yet"); |
sqfrsf(f,dcp); |
} |
} |
|
|
void sfusqfr(f,dcp) |
void sfusqfr(P f,DCP *dcp) |
P f; |
|
DCP *dcp; |
|
{ |
{ |
DCP dc,dct; |
DCP dc,dct; |
struct oDUM *udc; |
struct oDUM *udc; |
|
|
*dcp = dct; |
*dcp = dct; |
} |
} |
|
|
void sfbsqfr(f,x,y,dcp) |
#if 0 |
P f; |
void sfbsqfrmain(P f,V x,V y,DCP *dcp) |
V x,y; |
|
DCP *dcp; |
|
{ |
{ |
|
/* XXX*/ |
|
} |
|
|
|
/* f is bivariate */ |
|
|
|
void sfbsqfr(P f,V x,V y,DCP *dcp) |
|
{ |
P t,rf,cx,cy; |
P t,rf,cx,cy; |
VL vl,rvl; |
VL vl,rvl; |
DCP dcx,dcy; |
DCP dcx,dcy,dct,dc; |
struct oVL vl0,vl1; |
struct oVL vl0,vl1; |
|
|
/* vl = [x,y] */ |
/* cy(y) = cont(f,x), f /= cy */ |
vl0.v = x; vl0.next = &vl1; vl1.v = y; vl1.next = 0; vl = &vl0; |
|
/* cy(y) = cont(f,x), f /= cx */ |
|
cont_pp_sfp(vl,f,&cy,&t); f = t; |
cont_pp_sfp(vl,f,&cy,&t); f = t; |
/* rvl = [y,x] */ |
/* rvl = [y,x] */ |
reordvar(vl,y,&rvl); reorderp(rvl,vl,f,&rf); |
reordvar(vl,y,&rvl); reorderp(rvl,vl,f,&rf); |
|
|
reorderp(vl,rvl,rf,&f); |
reorderp(vl,rvl,rf,&f); |
|
|
/* f -> cx*cy*f */ |
/* f -> cx*cy*f */ |
sfusqfr(cx,&dcx); |
sfsqfr(cx,&dcx); dcx = NEXT(dcx); |
sfusqfr(cy,&dcy); |
sfsqfr(cy,&dcy); dcy = NEXT(dcy); |
|
if ( dcx ) { |
|
for ( dct = dcx; NEXT(dct); dct = NEXT(dct) ); |
|
NEXT(dct) = dcy; |
|
} else |
|
dcx = dcy; |
|
if ( OID(f) == O_N ) |
|
*dcp = dcx; |
|
else { |
|
/* f must be bivariate */ |
|
sfbsqfrmain(f,x,y,&dc); |
|
if ( dcx ) { |
|
for ( dct = dcx; NEXT(dct); dct = NEXT(dct) ); |
|
NEXT(dct) = dc; |
|
} else |
|
dcx = dc; |
|
*dcp = dcx; |
|
} |
} |
} |
|
#endif |
|
|
void sfdtest(P,ML,V,V,DCP *); |
void sfdtest(P,ML,V,V,DCP *); |
|
|
void sfbfctr(f,x,y,dcp) |
/* if degbound >= 0 find factor s.t. deg_x(factor) <= degbound */ |
P f; |
|
V x,y; |
void sfbfctr(P f,V x,V y,int degbound,DCP *dcp) |
DCP *dcp; |
|
{ |
{ |
ML list; |
ML list; |
P sf; |
P sf; |
|
|
int dx,dy; |
int dx,dy; |
|
|
/* sf(x) = f(x+ev) = list->c[0]*list->c[1]*... */ |
/* sf(x) = f(x+ev) = list->c[0]*list->c[1]*... */ |
sfhensel(5,f,x,&ev,&sf,&list); |
sfhensel(5,f,x,degbound,&ev,&sf,&list); |
if ( list->n == 0 ) |
if ( list->n == 0 ) |
error("sfbfctr : short of evaluation points"); |
error("sfbfctr : short of evaluation points"); |
else if ( list->n == 1 ) { |
else if ( list->n == 1 ) { |
|
|
*dcp = dc; |
*dcp = dc; |
} |
} |
|
|
|
/* returns shifted f, shifted factors and the eval pt */ |
|
|
|
void sfbfctr_shift(P f,V x,V y,int degbound,GFS *evp,P *sfp,DCP *dcp) |
|
{ |
|
ML list; |
|
P sf; |
|
GFS ev; |
|
DCP dc,dct; |
|
int dx,dy; |
|
|
|
/* sf(x) = f(x+ev) = list->c[0]*list->c[1]*... */ |
|
sfhensel(5,f,x,degbound,&ev,&sf,&list); |
|
if ( list->n == 0 ) |
|
error("sfbfctr_shift : short of evaluation points"); |
|
else if ( list->n == 1 ) { |
|
/* f is irreducible */ |
|
NEWDC(dc); DEG(dc) = ONE; COEF(dc) = f; NEXT(dc) = 0; |
|
*evp = 0; |
|
*sfp = f; |
|
*dcp = dc; |
|
} else { |
|
sfdtest(sf,list,x,y,dcp); |
|
*evp = ev; |
|
*sfp = sf; |
|
} |
|
} |
|
|
/* f = f(x,y) = list->c[0]*list->c[1]*... mod y^(list->bound+1) */ |
/* f = f(x,y) = list->c[0]*list->c[1]*... mod y^(list->bound+1) */ |
|
|
void sfdtest(f,list,x,y,dcp) |
void sfdtest(P f,ML list,V x,V y,DCP *dcp) |
P f; |
|
ML list; |
|
V x,y; |
|
DCP *dcp; |
|
{ |
{ |
int np,dx,dy; |
int np,dx,dy; |
int i,j,k; |
int i,j,k,bound; |
int *win; |
int *win; |
P g,lcg,factor,cofactor,lcyx; |
P g,lcg,factor,cofactor,lcyx; |
P t,csum; |
P csum; |
DCP dcf,dcf0,dc; |
DCP dcf,dcf0,dc; |
BM *c; |
BM *c; |
BM lcy; |
BM lcy; |
UM lcg0; |
UM lcg0,lcy0,w; |
|
UM *d1c; |
ML wlist; |
ML wlist; |
struct oVL vl1,vl0; |
struct oVL vl1,vl0; |
VL vl; |
VL vl; |
int z; |
int z,dt,dtok; |
|
|
/* vl = [x,y] */ |
/* vl = [x,y] */ |
vl0.v = x; vl0.next = &vl1; vl1.v = y; vl1.next = 0; vl = &vl0; |
vl0.v = x; vl0.next = &vl1; vl1.v = y; vl1.next = 0; vl = &vl0; |
|
|
win = W_ALLOC(np+1); |
win = W_ALLOC(np+1); |
wlist = W_MLALLOC(np); |
wlist = W_MLALLOC(np); |
wlist->n = list->n; |
wlist->n = list->n; |
wlist->bound = dy; |
bound = wlist->bound = list->bound; |
c = (BM *)COEF(wlist); |
c = (BM *)COEF(wlist); |
bcopy((char *)COEF(list),(char *)c,(int)(sizeof(BM)*np)); |
bcopy((char *)COEF(list),(char *)c,(int)(sizeof(BM)*np)); |
|
|
|
|
NEWP(lcyx); VR(lcyx) = x; DC(lcyx) = dc; |
NEWP(lcyx); VR(lcyx) = x; DC(lcyx) = dc; |
ptosfbm(dy,lcyx,lcy); |
ptosfbm(dy,lcyx,lcy); |
|
|
|
/* initialize lcy0 by LC(f) */ |
|
lcy0 = W_UMALLOC(bound); |
|
ptosfum(COEF(DC(g)),lcy0); |
|
|
|
/* ((d-1 coefs)*lcy0 */ |
|
d1c = (UM *)W_ALLOC(np*sizeof(UM)); |
|
w = W_UMALLOC(2*bound); |
|
for ( i = 1; i < np; i++ ) { |
|
extractcoefbm(c[i],degbm(c[i])-1,w); |
|
d1c[i] = W_UMALLOC(2*bound); |
|
mulsfum(w,lcy0,d1c[i]); |
|
/* d1c[i] = d1c[i] mod y^(bound+1) */ |
|
if ( DEG(d1c[i]) > bound ) { |
|
for ( j = DEG(d1c[i]); j > bound; j-- ) |
|
COEF(d1c[i])[j] = 0; |
|
degum(d1c[i],bound); |
|
} |
|
} |
|
|
fprintf(stderr,"np = %d\n",np); |
fprintf(stderr,"np = %d\n",np); |
|
dtok = 0; |
for ( g = f, k = 1, dcf = dcf0 = 0, win[0] = 1, --np, z = 0; ; z++ ) { |
for ( g = f, k = 1, dcf = dcf0 = 0, win[0] = 1, --np, z = 0; ; z++ ) { |
if ( !(z % 1000) ) fprintf(stderr,"."); |
if ( !(z % 1000) ) fprintf(stderr,"."); |
if ( sfdtestmain(vl,lcg,lcg0,lcy,csum,wlist,k,win,&factor,&cofactor) ) { |
dt = sfdegtest(dy,bound,d1c,k,win); |
|
if ( dt ) |
|
dtok++; |
|
if ( dt && sfdtestmain(vl,lcg,lcg0,lcy,csum,wlist, |
|
k,win,&factor,&cofactor) ) { |
NEXTDC(dcf0,dcf); DEG(dcf) = ONE; COEF(dcf) = factor; |
NEXTDC(dcf0,dcf); DEG(dcf) = ONE; COEF(dcf) = factor; |
g = cofactor; |
g = cofactor; |
|
|
|
|
/* update csum */ |
/* update csum */ |
sfcsump(vl,lcg,&csum); |
sfcsump(vl,lcg,&csum); |
|
|
|
/* update dy */ |
|
dy = getdeg(y,g); |
|
|
/* update lcy */ |
/* update lcy */ |
clearbm(0,lcy); |
clearbm(0,lcy); |
COEF(dc) = COEF(DC(g)); |
COEF(dc) = COEF(DC(g)); |
|
|
else |
else |
for ( i = 1; i < k; i++ ) |
for ( i = 1; i < k; i++ ) |
win[i] = win[0] + i; |
win[i] = win[0] + i; |
|
|
|
|
|
/* update lcy0 */ |
|
ptosfum(COEF(DC(g)),lcy0); |
|
|
|
/* update d-1 coeffs */ |
|
for ( i = 1; i <= np; i++ ) { |
|
extractcoefbm(c[i],degbm(c[i])-1,w); |
|
mulsfum(w,lcy0,d1c[i]); |
|
/* d1c[i] = d1c[1] mod y^(bound+1) */ |
|
if ( DEG(d1c[i]) > bound ) { |
|
for ( j = DEG(d1c[i]); j > bound; j-- ) |
|
COEF(d1c[i])[j] = 0; |
|
degum(d1c[i],bound); |
|
} |
|
} |
} else if ( !ncombi(1,np,k,win) ) |
} else if ( !ncombi(1,np,k,win) ) |
if ( k == np ) |
if ( k == np ) |
break; |
break; |
|
|
for ( i = 0, ++k; i < k; i++ ) |
for ( i = 0, ++k; i < k; i++ ) |
win[i] = i + 1; |
win[i] = i + 1; |
} |
} |
fprintf(stderr,"\n"); |
fprintf(stderr,"total %d, omitted by degtest %d\n",z,z-dtok); |
NEXTDC(dcf0,dcf); COEF(dcf) = g; |
NEXTDC(dcf0,dcf); COEF(dcf) = g; |
DEG(dcf) = ONE; NEXT(dcf) = 0; *dcp = dcf0; |
DEG(dcf) = ONE; NEXT(dcf) = 0; *dcp = dcf0; |
} |
} |
|
|
|
void extractcoefbm(BM f,int dx,UM r) |
|
{ |
|
int j; |
|
UM fj; |
|
|
|
for ( j = DEG(f); j >= 0; j-- ) { |
|
fj = COEF(f)[j]; |
|
if ( fj && DEG(fj) >= dx ) { |
|
COEF(r)[j] = COEF(fj)[dx]; |
|
} else |
|
COEF(r)[j] = 0; |
|
} |
|
degum(r,DEG(f)); |
|
} |
|
|
|
/* deg_y(prod mod y^(bound+1)) <= dy ? */ |
|
|
|
int sfdegtest(int dy,int bound,UM *d1c,int k,int *in) |
|
{ |
|
int i,j; |
|
UM w,w1,wt; |
|
BM t; |
|
|
|
w = W_UMALLOC(bound); |
|
w1 = W_UMALLOC(bound); |
|
clearum(w,bound); |
|
for ( i = 0; i < k; i++ ) { |
|
addsfum(w,d1c[in[i]],w1); wt = w; w = w1; w1 = wt; |
|
} |
|
return DEG(w) <= dy ? 1 : 0; |
|
} |
|
|
/* lcy = LC(g), lcg = lcy*g, lcg0 = const part of lcg */ |
/* lcy = LC(g), lcg = lcy*g, lcg0 = const part of lcg */ |
int sfdtestmain(vl,lcg,lcg0,lcy,csum,list,k,in,fp,cofp) |
|
VL vl; |
int sfdtestmain(VL vl,P lcg,UM lcg0,BM lcy,P csum,ML list, |
P lcg; |
int k,int *in,P *fp,P *cofp) |
UM lcg0; |
|
BM lcy; |
|
P csum; |
|
ML list; |
|
int k; |
|
int *in; |
|
P *fp,*cofp; |
|
{ |
{ |
P fmul,csumg,q,cont; |
P fmul,csumg,q,cont; |
V x,y; |
V x,y; |
|
|
return 0; |
return 0; |
} |
} |
|
|
void const_term(f,c) |
void const_term(P f,UM c) |
P f; |
|
UM c; |
|
{ |
{ |
DCP dc; |
DCP dc; |
|
|
|
|
DEG(c) = -1; |
DEG(c) = -1; |
} |
} |
|
|
void const_term_sfbm(f,c) |
void const_term_sfbm(BM f,UM c) |
BM f; |
|
UM c; |
|
{ |
{ |
int i,dy; |
int i,dy; |
|
|
|
|
|
|
/* lcy*(product of const part) | lcg0 ? */ |
/* lcy*(product of const part) | lcg0 ? */ |
|
|
int sfctest(lcg0,lcy,list,k,in) |
int sfctest(UM lcg0,BM lcy,ML list,int k,int *in) |
UM lcg0; |
|
BM lcy; |
|
ML list; |
|
int k; |
|
int *in; |
|
{ |
{ |
DCP dc; |
|
int dy,i,dr; |
int dy,i,dr; |
UM t,s,u,w; |
UM t,s,u,w; |
BM *l; |
BM *l; |
|
|
|
|
/* main var of f is x */ |
/* main var of f is x */ |
|
|
void mulsfbmarray(dx,lcy,list,k,in,x,y,g) |
void mulsfbmarray(int dx,BM lcy,ML list,int k,int *in,V x,V y,P *g) |
int dx; |
|
BM lcy; |
|
ML list; |
|
int k; |
|
int *in; |
|
V x,y; |
|
P *g; |
|
{ |
{ |
int dy,i; |
int dy,i; |
BM wb0,wb1,t,lcbm; |
BM wb0,wb1,t; |
BM *l; |
BM *l; |
|
|
dy = list->bound; |
dy = list->bound; |
|
|
sfbmtop(wb0,x,y,g); |
sfbmtop(wb0,x,y,g); |
} |
} |
|
|
void sfcsump(vl,f,s) |
void sfcsump(VL vl,P f,P *s) |
VL vl; |
|
P f; |
|
P *s; |
|
{ |
{ |
P t,u; |
P t,u; |
DCP dc; |
DCP dc; |
|
|
|
|
/* *fp = primitive part of f w.r.t. x */ |
/* *fp = primitive part of f w.r.t. x */ |
|
|
void cont_pp_sfp(vl,f,cp,fp) |
void cont_pp_sfp(VL vl,P f,P *cp,P *fp) |
VL vl; |
|
P f; |
|
P *cp,*fp; |
|
{ |
{ |
V x,y; |
V x,y; |
int d; |
int d; |
|
|
y = vl->next->v; |
y = vl->next->v; |
d = getdeg(y,f); |
d = getdeg(y,f); |
if ( d == 0 ) { |
if ( d == 0 ) { |
MKGFS(0,g); |
itogfs(1,&g); |
*cp = (P)g; |
*cp = (P)g; |
*fp = f; /* XXX */ |
*fp = f; /* XXX */ |
} else { |
} else { |
|
|
} |
} |
} |
} |
|
|
int divtp_by_sfbm(vl,f,g,qp) |
int divtp_by_sfbm(VL vl,P f,P g,P *qp) |
VL vl; |
|
P f,g; |
|
P *qp; |
|
{ |
{ |
V x,y; |
V x,y; |
int fx,fy,gx,gy; |
int fx,fy,gx,gy; |
|
|
if ( j >= 0 ) |
if ( j >= 0 ) |
return 0; |
return 0; |
sfbmtop(ql,x,y,qp); |
sfbmtop(ql,x,y,qp); |
|
return 1; |
} |
} |
|
|
/* XXX generate an irreducible poly of degree n */ |
/* XXX generate an irreducible poly of degree n */ |
|
|
extern int current_gfs_q1; |
extern int current_gfs_q1; |
|
extern int *current_gfs_ntoi; |
|
|
void generate_defpoly_sfum(n,dp) |
void generate_defpoly_sfum(int n,UM *dp) |
int n; |
|
UM *dp; |
|
{ |
{ |
UM r,dr,t,g; |
UM r,dr,t,g; |
UM *f; |
UM *f; |
|
|
for ( j = 0; j < i; j++ ) |
for ( j = 0; j < i; j++ ) |
w[j] = 0; |
w[j] = 0; |
w[i]++; |
w[i]++; |
for ( i = 0; i < n; i++ ) |
if ( !current_gfs_ntoi ) |
c[i] = w[i]?FTOIF(w[i]-1):0; |
for ( i = 0; i < n; i++ ) |
|
c[i] = w[i]?FTOIF(w[i]):0; |
|
else |
|
for ( i = 0; i < n; i++ ) |
|
c[i] = w[i]?FTOIF(w[i]-1):0; |
if ( !c[0] ) |
if ( !c[0] ) |
continue; |
continue; |
diffsfum(r,dr); cpyum(r,t); gcdsfum(t,dr,g); |
diffsfum(r,dr); cpyum(r,t); gcdsfum(t,dr,g); |