=================================================================== RCS file: /home/cvs/OpenXM_contrib2/asir2000/engine/nd.c,v retrieving revision 1.119 retrieving revision 1.122 diff -u -p -r1.119 -r1.122 --- OpenXM_contrib2/asir2000/engine/nd.c 2004/12/06 09:29:34 1.119 +++ OpenXM_contrib2/asir2000/engine/nd.c 2004/12/09 13:20:33 1.122 @@ -1,4 +1,4 @@ -/* $OpenXM: OpenXM_contrib2/asir2000/engine/nd.c,v 1.118 2004/12/04 09:39:27 noro Exp $ */ +/* $OpenXM: OpenXM_contrib2/asir2000/engine/nd.c,v 1.121 2004/12/09 08:56:43 noro Exp $ */ #include "nd.h" @@ -8,6 +8,7 @@ NM _nm_free_list; ND _nd_free_list; ND_pairs _ndp_free_list; +static int nd_ntrans; static int nd_nalg; #if 0 static int ndv_alloc; @@ -1536,7 +1537,7 @@ ND normalize_pbucket(int mod,PGeoBucket g) return r; } -void do_diagonalize(int sugar,int m) +int do_diagonalize(int sugar,int m) { int i,nh,stat; NODE r,g,t; @@ -1552,7 +1553,8 @@ void do_diagonalize(int sugar,int m) nfv = nd_ps[i]; s = ndvtond(m,nfv); s = nd_separate_head(s,&head); - nd_nf(m,s,nd_ps,1,&dn,&nf); + stat = nd_nf(m,s,nd_ps,1,&dn,&nf); + if ( !stat ) return 0; if ( !m ) { NTOQ(NM(dn.z),SGN(dn.z),num); mulq(HCQ(head),num,&q); HCQ(head) = q; @@ -1573,6 +1575,7 @@ void do_diagonalize(int sugar,int m) } else nd_ps[i] = nfv; } + return 1; } /* return value = 0 => input is not a GB */ @@ -1598,8 +1601,14 @@ NODE nd_gb(int m,int ishomo,int checkonly) again: l = nd_minp(d,&d); if ( SG(l) != sugar ) { - if ( ishomo ) do_diagonalize(sugar,m); - + if ( ishomo ) { + stat = do_diagonalize(sugar,m); + if ( !stat ) { + NEXT(l) = d; d = l; + d = nd_reconstruct(0,d); + goto again; + } + } sugar = SG(l); if ( DP_Print ) fprintf(asir_out,"%d",sugar); } @@ -1645,7 +1654,7 @@ again: return g; } -void do_diagonalize_trace(int sugar,int m) +int do_diagonalize_trace(int sugar,int m) { int i,nh,stat; NODE r,g,t; @@ -1658,7 +1667,8 @@ void do_diagonalize_trace(int sugar,int m) /* for nd_ps */ s = ndvtond(m,nd_ps[i]); s = nd_separate_head(s,&head); - nd_nf_pbucket(m,s,nd_ps,1,&nf); + stat = nd_nf_pbucket(m,s,nd_ps,1,&nf); + if ( !stat ) return 0; nf = nd_add(m,head,nf); ndv_free(nd_ps[i]); nd_ps[i] = ndtondv(m,nf); @@ -1671,7 +1681,8 @@ void do_diagonalize_trace(int sugar,int m) nfv = nd_ps_trace[i]; s = ndvtond(0,nfv); s = nd_separate_head(s,&head); - nd_nf(0,s,nd_ps_trace,1,&dn,&nf); + stat = nd_nf(0,s,nd_ps_trace,1,&dn,&nf); + if ( !stat ) return 0; NTOQ(NM(dn.z),SGN(dn.z),num); mulq(HCQ(head),num,&q); HCQ(head) = q; if ( DN(dn.z) ) { @@ -1690,6 +1701,7 @@ void do_diagonalize_trace(int sugar,int m) } else nd_ps_trace[i] = nfv; } + return 1; } static struct oEGT eg_invdalg; @@ -1720,7 +1732,14 @@ NODE nd_gb_trace(int m,int ishomo) again: l = nd_minp(d,&d); if ( SG(l) != sugar ) { - if ( ishomo ) do_diagonalize_trace(sugar,m); + if ( ishomo ) { + stat = do_diagonalize_trace(sugar,m); + if ( !stat ) { + NEXT(l) = d; d = l; + d = nd_reconstruct(1,d); + goto again; + } + } sugar = SG(l); if ( DP_Print ) fprintf(asir_out,"%d",sugar); } @@ -2214,6 +2233,108 @@ void ndv_setup(int mod,int trace,NODE f) struct order_spec *append_block(struct order_spec *spec, int nv,int nalg,int ord); +extern VECT current_dl_weight_vector_obj; +static VECT prev_weight_vector_obj; + +void preprocess_algcoef(VL vv,VL av,struct order_spec *ord,LIST f, + struct order_spec **ord1p,LIST *f1p,NODE *alistp) +{ + NODE alist,t,s,r0,r,arg; + VL tv; + P poly; + DP d; + Alg alpha,dp; + DAlg inv,da,hc; + MP m; + int i,nvar,nalg,n; + NumberField nf; + LIST f1,f2; + struct order_spec *current_spec; + VECT obj,obj0; + Obj tmp; + + for ( nvar = 0, tv = vv; tv; tv = NEXT(tv), nvar++); + for ( nalg = 0, tv = av; tv; tv = NEXT(tv), nalg++); + + for ( alist = 0, tv = av; tv; tv = NEXT(tv) ) { + NEXTNODE(alist,t); MKV(tv->v,poly); + MKAlg(poly,alpha); BDY(t) = (pointer)alpha; + tv->v = tv->v->priv; + } + NEXT(t) = 0; + + /* simplification, makeing polynomials monic */ + setfield_dalg(alist); + obj_algtodalg(f,&f1); + for ( t = BDY(f); t; t = NEXT(t) ) { + initd(ord); ptod(vv,vv,(P)BDY(t),&d); + hc = (DAlg)BDY(d)->c; + if ( NID(hc) == N_DA ) { + invdalg(hc,&inv); + for ( m = BDY(d); m; m = NEXT(m) ) { + muldalg(inv,(DAlg)m->c,&da); m->c = (P)da; + } + } + initd(ord); dtop(vv,vv,d,&poly); BDY(f) = (pointer)poly; + } + obj_dalgtoalg(f1,&f); + + /* append alg vars to the var list */ + for ( tv = vv; NEXT(tv); tv = NEXT(tv) ); + NEXT(tv) = av; + + /* append a block to ord */ + *ord1p = append_block(ord,nvar,nalg,2); + + /* create generator list */ + nf = get_numberfield(); + for ( i = nalg-1, t = BDY(f); i >= 0; i-- ) { + MKAlg(nf->defpoly[i],dp); + MKNODE(s,dp,t); t = s; + } + MKLIST(f1,t); + *alistp = alist; + algobjtorat(f1,f1p); + + /* creating a new weight vector */ + prev_weight_vector_obj = obj0 = current_dl_weight_vector_obj; + n = nvar+nalg+1; + MKVECT(obj,n); + if ( obj0 && obj0->len == nvar ) + for ( i = 0; i < nvar; i++ ) BDY(obj)[i] = BDY(obj0)[i]; + else + for ( i = 0; i < nvar; i++ ) BDY(obj)[i] = (pointer)ONE; + for ( i = 0; i < nalg; i++ ) BDY(obj)[i+nvar] = 0; + BDY(obj)[n-1] = (pointer)ONE; + arg = mknode(1,obj); + Pdp_set_weight(arg,&tmp); +} + +NODE postprocess_algcoef(VL av,NODE alist,NODE r) +{ + NODE s,t,u0,u; + P p; + VL tv; + Obj obj,tmp; + NODE arg; + + u0 = 0; + for ( t = r; t; t = NEXT(t) ) { + p = (P)BDY(t); + for ( tv = av, s = alist; tv; tv = NEXT(tv), s = NEXT(s) ) { + substr(CO,0,(Obj)p,tv->v,(Obj)BDY(s),&obj); p = (P)obj; + } + if ( OID(p) == O_P || (OID(p) == O_N && NID((Num)p) != N_A) ) { + NEXTNODE(u0,u); + BDY(u) = (pointer)p; + } + } + arg = mknode(1,prev_weight_vector_obj); + Pdp_set_weight(arg,&tmp); + + return u0; +} + void nd_gr(LIST f,LIST v,int m,int f4,struct order_spec *ord,LIST *rp) { VL tv,fv,vv,vc,av; @@ -2226,6 +2347,7 @@ void nd_gr(LIST f,LIST v,int m,int f4,struct order_spe LIST f1,f2; Obj obj; NumberField nf; + struct order_spec *ord1; if ( !m && Demand ) nd_demand = 1; else nd_demand = 0; @@ -2247,28 +2369,15 @@ void nd_gr(LIST f,LIST v,int m,int f4,struct order_spe if ( !m ) { get_algtree((Obj)f,&av); for ( nalg = 0, tv = av; tv; tv = NEXT(tv), nalg++ ); + nd_ntrans = nvar; nd_nalg = nalg; /* #i -> t#i */ if ( nalg ) { - for ( alist = 0, tv = av; tv; tv = NEXT(tv) ) { - NEXTNODE(alist,t); MKV(tv->v,p); - MKAlg(p,alpha); BDY(t) = (pointer)alpha; - tv->v = tv->v->priv; - } - NEXT(t) = 0; - for ( tv = vv; NEXT(tv); tv = NEXT(tv) ); - NEXT(tv) = av; - ord = append_block(ord,nvar,nalg,2); - nvar += nalg; - setfield_dalg(alist); - nf = get_numberfield(); - for ( i = nalg-1, t = BDY(f); i >= 0; i-- ) { - MKAlg(nf->defpoly[i],dp); - MKNODE(s,dp,t); t = s; - } - MKLIST(f1,t); - algobjtorat(f1,&f2); f = f2; + preprocess_algcoef(vv,av,ord,f,&ord1,&f1,&alist); + ord = ord1; + f = f1; } + nvar += nalg; } nd_init_ord(ord); for ( t = BDY(f), max = 0; t; t = NEXT(t) ) @@ -2294,15 +2403,10 @@ void nd_gr(LIST f,LIST v,int m,int f4,struct order_spe for ( r0 = 0, t = x; t; t = NEXT(t) ) { NEXTNODE(r0,r); BDY(r) = ndvtop(m,CO,vv,BDY(t)); - if ( nalg ) { - p = BDY(r); - for ( tv = av, s = alist; tv; tv = NEXT(tv), s = NEXT(s) ) { - substr(CO,0,(Obj)p,tv->v,(Obj)BDY(s),&obj); p = (P)obj; - } - BDY(r) = p; - } } if ( r0 ) NEXT(r) = 0; + if ( nalg ) + r0 = postprocess_algcoef(av,alist,r0); MKLIST(*rp,r0); #if 0 fprintf(asir_out,"ndv_alloc=%d\n",ndv_alloc); @@ -2311,7 +2415,6 @@ void nd_gr(LIST f,LIST v,int m,int f4,struct order_spe void nd_gr_trace(LIST f,LIST v,int trace,int homo,struct order_spec *ord,LIST *rp) { - struct order_spec *ord1; VL tv,fv,vv,vc,av; NODE fd,fd0,in0,in,r,r0,t,s,cand,alist; int m,nocheck,nvar,mindex,e,max; @@ -2325,6 +2428,7 @@ void nd_gr_trace(LIST f,LIST v,int trace,int homo,stru LIST f1,f2; Obj obj; NumberField nf; + struct order_spec *ord1; get_vars((Obj)f,&fv); pltovl(v,&vv); for ( nvar = 0, tv = vv; tv; tv = NEXT(tv), nvar++ ); @@ -2339,28 +2443,15 @@ void nd_gr_trace(LIST f,LIST v,int trace,int homo,stru get_algtree((Obj)f,&av); for ( nalg = 0, tv = av; tv; tv = NEXT(tv), nalg++ ); + nd_ntrans = nvar; nd_nalg = nalg; /* #i -> t#i */ if ( nalg ) { - for ( alist = 0, tv = av; tv; tv = NEXT(tv) ) { - NEXTNODE(alist,t); MKV(tv->v,poly); - MKAlg(poly,alpha); BDY(t) = (pointer)alpha; - tv->v = tv->v->priv; - } - NEXT(t) = 0; - for ( tv = vv; NEXT(tv); tv = NEXT(tv) ); - NEXT(tv) = av; - ord = append_block(ord,nvar,nalg,2); - nvar += nalg; - setfield_dalg(alist); - nf = get_numberfield(); - for ( i = nalg-1, t = BDY(f); i >= 0; i-- ) { - MKAlg(nf->defpoly[i],dp); - MKNODE(s,dp,t); t = s; - } - MKLIST(f1,t); - algobjtorat(f1,&f2); f = f2; + preprocess_algcoef(vv,av,ord,f,&ord1,&f1,&alist); + ord = ord1; + f = f1; } + nvar += nalg; nocheck = 0; mindex = 0; @@ -2395,7 +2486,7 @@ void nd_gr_trace(LIST f,LIST v,int trace,int homo,stru if ( in0 ) NEXT(in) = 0; if ( fd0 ) NEXT(fd) = 0; if ( !ishomo && homo ) { - for ( t = in0, wmax = 0; t; t = NEXT(t) ) { + for ( t = in0, wmax = max; t; t = NEXT(t) ) { c = (NDV)BDY(t); len = LEN(c); for ( a = BDY(c), i = 0; i < len; i++, NMV_ADV(a) ) wmax = MAX(TD(DL(a)),wmax); @@ -2451,16 +2542,10 @@ void nd_gr_trace(LIST f,LIST v,int trace,int homo,stru /* dp->p */ nd_bpe = cbpe; nd_setup_parameters(nd_nvar,0); - for ( r = cand; r; r = NEXT(r) ) { + for ( r = cand; r; r = NEXT(r) ) BDY(r) = (pointer)ndvtop(0,CO,vv,BDY(r)); - if ( nalg ) { - poly = BDY(r); - for ( tv = av, s = alist; tv; tv = NEXT(tv), s = NEXT(s) ) { - substr(CO,0,(Obj)poly,tv->v,(Obj)BDY(s),&obj); poly = (P)obj; - } - BDY(r) = poly; - } - } + if ( nalg ) + cand = postprocess_algcoef(av,alist,cand); MKLIST(*rp,cand); } @@ -5209,24 +5294,31 @@ UINT *nd_det_compute_bound(NDV **dm,int n,int j) DL nd_separate_d(UINT *d,UINT *trans) { - int n,ntrans,td,i,e; + int n,td,i,e,j; DL a; - n = nd_nvar; ntrans = n-nd_nalg; ndl_zero(trans); td = 0; - for ( i = 0; i < ntrans; i++ ) { + for ( i = 0; i < nd_ntrans; i++ ) { e = GET_EXP(d,i); PUT_EXP(trans,i,e); td += MUL_WEIGHT(e,i); } + if ( nd_ntrans+nd_nalg < nd_nvar ) { + /* homogenized */ + i = nd_nvar-1; + e = GET_EXP(d,i); + PUT_EXP(trans,i,e); + td += MUL_WEIGHT(e,i); + } TD(trans) = td; if ( nd_blockmask) ndl_weight_mask(trans); NEWDL(a,nd_nalg); td = 0; - for ( ; i < n; i++ ) { - e = GET_EXP(d,i); - a->d[i-ntrans] = e; + for ( i = 0; i < nd_nalg; i++ ) { + j = nd_ntrans+i; + e = GET_EXP(d,j); + a->d[i] = e; td += e; } a->td = td; @@ -5254,7 +5346,6 @@ int nd_monic(int mod,ND *p) if ( !(nf = get_numberfield()) ) error("nd_monic : current_numberfield is not set"); - n = nd_nvar; ntrans = n-nd_nalg; /* Q coef -> DAlg coef */ NEWNM(ma0); ma = ma0; m = BDY(*p); @@ -5308,10 +5399,10 @@ int nd_monic(int mod,ND *p) dl = mp->dl; td = TD(DL(m)); ndl_copy(DL(m),DL(mr)); - for ( i = ntrans; i < n; i++ ) { - e = dl->d[i-ntrans]; - PUT_EXP(DL(mr),i,e); - td += MUL_WEIGHT(e,i); + for ( i = 0; i < nd_nalg; i++ ) { + e = dl->d[i]; + PUT_EXP(DL(mr),i+nd_ntrans,e); + td += MUL_WEIGHT(e,i+nd_ntrans); } TD(DL(mr)) = td; if ( nd_blockmask) ndl_weight_mask(DL(mr));