=================================================================== RCS file: /home/cvs/OpenXM_contrib2/asir2000/engine/dist.c,v retrieving revision 1.12 retrieving revision 1.16 diff -u -p -r1.12 -r1.16 --- OpenXM_contrib2/asir2000/engine/dist.c 2000/12/11 02:00:41 1.12 +++ OpenXM_contrib2/asir2000/engine/dist.c 2001/05/02 09:03:53 1.16 @@ -45,7 +45,7 @@ * DEVELOPER SHALL HAVE NO LIABILITY IN CONNECTION WITH THE USE, * PERFORMANCE OR NON-PERFORMANCE OF THE SOFTWARE. * - * $OpenXM: OpenXM_contrib2/asir2000/engine/dist.c,v 1.11 2000/12/05 06:59:16 noro Exp $ + * $OpenXM: OpenXM_contrib2/asir2000/engine/dist.c,v 1.15 2001/03/29 09:49:57 noro Exp $ */ #include "ca.h" @@ -60,6 +60,7 @@ #define ORD_BLEXREV 8 #define ORD_ELIM 9 #define ORD_WEYL_ELIM 10 +#define ORD_HOMO_WW_DRL 11 int (*cmpdl)()=cmpdl_revgradlex; int (*primitive_cmpdl[3])() = {cmpdl_revgradlex,cmpdl_gradlex,cmpdl_lex}; @@ -73,6 +74,7 @@ void comm_muld_tab(VL,int,struct cdl *,int,struct cdl void mkwc(int,int,Q *); int cmpdl_weyl_elim(); +int cmpdl_homo_ww_drl(); int do_weyl; @@ -104,7 +106,10 @@ P f; if ( !f ) return 0; else if ( NUM(f) ) - return (NID((Num)f) == N_LM || NID((Num)f) == N_GF2N) ? 1 : 0; + return (NID((Num)f) == N_LM + || NID((Num)f) == N_GF2N + || NID((Num)f) == N_GFPN + || NID((Num)f) == N_GFS) ? 1 : 0; else { for ( dc = DC(f); dc; dc = NEXT(dc) ) if ( has_fcoef_p(COEF(dc)) ) @@ -146,6 +151,8 @@ struct order_spec *spec; cmpdl = cmpdl_elim; break; case ORD_WEYL_ELIM: cmpdl = cmpdl_weyl_elim; break; + case ORD_HOMO_WW_DRL: + cmpdl = cmpdl_homo_ww_drl; break; case ORD_LEX: default: cmpdl = cmpdl_lex; break; } @@ -160,12 +167,13 @@ P p; DP *pr; { int isconst = 0; - int n,i; + int n,i,j,k; VL tvl; V v; DL d; MP m; DCP dc; + DCP *w; DP r,s,t,u; P x,c; @@ -180,15 +188,25 @@ DP *pr; for ( i = 0, tvl = dvl, v = VR(p); tvl && tvl->v != v; tvl = NEXT(tvl), i++ ); if ( !tvl ) { - for ( dc = DC(p), s = 0, MKV(v,x); dc; dc = NEXT(dc) ) { - ptod(vl,dvl,COEF(dc),&t); pwrp(vl,x,DEG(dc),&c); + for ( dc = DC(p), k = 0; dc; dc = NEXT(dc), k++ ); + w = (DCP *)ALLOCA(k*sizeof(DCP)); + for ( dc = DC(p), j = 0; j < k; dc = NEXT(dc), j++ ) + w[j] = dc; + + for ( j = k-1, s = 0, MKV(v,x); j >= 0; j-- ) { + ptod(vl,dvl,COEF(w[j]),&t); pwrp(vl,x,DEG(w[j]),&c); muldc(vl,t,c,&r); addd(vl,r,s,&t); s = t; } *pr = s; } else { - for ( dc = DC(p), s = 0; dc; dc = NEXT(dc) ) { - ptod(vl,dvl,COEF(dc),&t); - NEWDL(d,n); d->td = QTOS(DEG(dc)); d->d[i] = d->td; + for ( dc = DC(p), k = 0; dc; dc = NEXT(dc), k++ ); + w = (DCP *)ALLOCA(k*sizeof(DCP)); + for ( dc = DC(p), j = 0; j < k; dc = NEXT(dc), j++ ) + w[j] = dc; + + for ( j = k-1, s = 0; j >= 0; j-- ) { + ptod(vl,dvl,COEF(w[j]),&t); + NEWDL(d,n); d->td = QTOS(DEG(w[j])); d->d[i] = d->td; NEWMP(m); m->dl = d; C(m) = (P)ONE; NEXT(m) = 0; MKDP(n,m,u); u->sugar = d->td; comm_muld(vl,t,u,&r); addd(vl,r,s,&t); s = t; } @@ -205,9 +223,10 @@ VL vl,dvl; DP p; P *pr; { - int n,i; + int n,i,j,k; DL d; MP m; + MP *a; P r,s,t,u,w; Q q; VL tvl; @@ -215,7 +234,13 @@ P *pr; if ( !p ) *pr = 0; else { - for ( n = p->nv, m = BDY(p), s = 0; m; m = NEXT(m) ) { + for ( k = 0, m = BDY(p); m; m = NEXT(m), k++ ); + a = (MP *)ALLOCA(k*sizeof(MP)); + for ( j = 0, m = BDY(p); j < k; m = NEXT(m), j++ ) + a[j] = m; + + for ( n = p->nv, j = k-1, s = 0; j >= 0; j-- ) { + m = a[j]; t = C(m); if ( NUM(t) && NID((Num)t) == N_M ) { mptop(t,&u); t = u; @@ -1092,6 +1117,49 @@ DL d1,d2; else if ( d1->td < d2->td ) return -1; else return -cmpdl_revlex(n,d1,d2); +} + +/* + a special ordering + 1. total order + 2. (-w,w) for the first 2*m variables + 3. DRL for the first 2*m variables +*/ + +extern int *current_weight_vector; + +int cmpdl_homo_ww_drl(n,d1,d2) +int n; +DL d1,d2; +{ + int e1,e2,m,i; + int *p1,*p2; + + if ( d1->td > d2->td ) + return 1; + else if ( d1->td < d2->td ) + return -1; + + m = n>>1; + for ( i = 0, e1 = e2 = 0; i < m; i++ ) { + e1 += current_weight_vector[i]*(d1->d[m+i] - d1->d[i]); + e2 += current_weight_vector[i]*(d2->d[m+i] - d2->d[i]); + } + if ( e1 > e2 ) + return 1; + else if ( e1 < e2 ) + return -1; + + e1 = d1->td - d1->d[n-1]; + e2 = d2->td - d2->d[n-1]; + if ( e1 > e2 ) + return 1; + else if ( e1 < e2 ) + return -1; + + for ( i= n - 1, p1 = d1->d+n-1, p2 = d2->d+n-1; + i >= 0 && *p1 == *p2; i--, p1--, p2-- ); + return i < 0 ? 0 : (*p1 < *p2 ? 1 : -1); } int cmpdl_order_pair(n,d1,d2)