=================================================================== RCS file: /home/cvs/OpenXM_contrib2/asir2000/engine/dist.c,v retrieving revision 1.26 retrieving revision 1.32 diff -u -p -r1.26 -r1.32 --- OpenXM_contrib2/asir2000/engine/dist.c 2003/08/22 08:14:45 1.26 +++ OpenXM_contrib2/asir2000/engine/dist.c 2004/06/15 16:14:50 1.32 @@ -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.25 2003/07/09 07:11:08 noro Exp $ + * $OpenXM: OpenXM_contrib2/asir2000/engine/dist.c,v 1.31 2004/05/14 06:02:54 noro Exp $ */ #include "ca.h" @@ -72,7 +72,8 @@ int (*primitive_cmpdl[3])() = {cmpdl_revgradlex,cmpdl_ int do_weyl; int dp_nelim,dp_fcoeffs; -struct order_spec dp_current_spec; +struct order_spec *dp_current_spec; +struct modorder_spec *dp_current_modspec; int *dp_dl_work; void comm_muld_trunc(VL vl,DP p1,DP p2,DL dl,DP *pr); @@ -80,6 +81,15 @@ void comm_quod(VL vl,DP p1,DP p2,DP *pr); void muldm_trunc(VL vl,DP p,MP m0,DL dl,DP *pr); void muldc_trunc(VL vl,DP p,P c,DL dl,DP *pr); +void order_init() +{ + struct order_spec *spec; + + create_order_spec(0,0,&spec); + initd(spec); + create_modorder_spec(0,0,&dp_current_modspec); +} + int has_sfcoef(DP f) { MP t; @@ -111,6 +121,10 @@ int has_sfcoef_p(P f) void initd(struct order_spec *spec) { switch ( spec->id ) { + case 3: + cmpdl = cmpdl_composite; + dp_dl_work = (int *)MALLOC_ATOMIC(spec->nv*sizeof(int)); + break; case 2: cmpdl = cmpdl_matrix; dp_dl_work = (int *)MALLOC_ATOMIC(spec->nv*sizeof(int)); @@ -151,7 +165,7 @@ void initd(struct order_spec *spec) } break; } - dp_current_spec = *spec; + dp_current_spec = spec; } void ptod(VL vl,VL dvl,P p,DP *pr) @@ -286,12 +300,23 @@ void addd(VL vl,DP p1,DP p2,DP *pr) int n; MP m1,m2,mr,mr0; P t; + DL d; if ( !p1 ) *pr = p2; else if ( !p2 ) *pr = p1; else { + if ( OID(p1) <= O_R ) { + n = NV(p2); NEWDL(d,n); + NEWMP(m1); m1->dl = d; C(m1) = (P)p1; NEXT(m1) = 0; + MKDP(n,m1,p1); (p1)->sugar = 0; + } + if ( OID(p2) <= O_R ) { + n = NV(p1); NEWDL(d,n); + NEWMP(m2); m2->dl = d; C(m2) = (P)p2; NEXT(m2) = 0; + MKDP(n,m2,p2); (p2)->sugar = 0; + } for ( n = NV(p1), m1 = BDY(p1), m2 = BDY(p2), mr0 = 0; m1 && m2; ) switch ( (*cmpdl)(n,m1->dl,m2->dl) ) { case 0: @@ -574,6 +599,8 @@ void chsgnd(DP p,DP *pr) if ( !p ) *pr = 0; + else if ( OID(p) <= O_R ) + chsgnr(p,pr); else { for ( mr0 = 0, m = BDY(p); m; m = NEXT(m) ) { NEXTMP(mr0,mr); chsgnp(C(m),&C(mr)); mr->dl = m->dl; @@ -1467,8 +1494,8 @@ int cmpdl_order_pair(int n,DL d1,DL d2) int len,head; struct order_pair *pair; - len = dp_current_spec.ord.block.length; - pair = dp_current_spec.ord.block.order_pair; + len = dp_current_spec->ord.block.length; + pair = dp_current_spec->ord.block.order_pair; head = 0; for ( i = 0, t1 = d1->d, t2 = d2->d; i < len; i++ ) { @@ -1517,6 +1544,75 @@ int cmpdl_order_pair(int n,DL d1,DL d2) return 0; } +int cmpdl_composite(int nv,DL d1,DL d2) +{ + int n,i,j,k,start,s,len; + int *dw; + struct sparse_weight *sw; + struct weight_or_block *worb; + int *w,*t1,*t2; + + n = dp_current_spec->ord.composite.length; + worb = dp_current_spec->ord.composite.w_or_b; + w = dp_dl_work; + for ( i = 0, t1 = d1->d, t2 = d2->d; i < nv; i++ ) + w[i] = t1[i]-t2[i]; + for ( i = 0; i < n; i++, worb++ ) { + len = worb->length; + switch ( worb->type ) { + case IS_DENSE_WEIGHT: + dw = worb->body.dense_weight; + for ( j = 0, s = 0; j < len; j++ ) + s += dw[j]*w[j]; + if ( s > 0 ) return 1; + else if ( s < 0 ) return -1; + break; + case IS_SPARSE_WEIGHT: + sw = worb->body.sparse_weight; + for ( j = 0, s = 0; j < len; j++ ) + s += sw[j].value*w[sw[j].pos]; + if ( s > 0 ) return 1; + else if ( s < 0 ) return -1; + break; + case IS_BLOCK: + start = worb->body.block.start; + switch ( worb->body.block.order ) { + case 0: + for ( j = 0, k = start, s = 0; j < len; j++, k++ ) { + s += MUL_WEIGHT(w[k],k); + } + if ( s > 0 ) return 1; + else if ( s < 0 ) return -1; + else { + for ( j = k-1; j >= start && w[j] == 0; j-- ); + if ( j >= start ) + return w[j] < 0 ? 1 : -1; + } + break; + case 1: + for ( j = 0, k = start, s = 0; j < len; j++, k++ ) { + s += MUL_WEIGHT(w[k],k); + } + if ( s > 0 ) return 1; + else if ( s < 0 ) return -1; + else { + for ( j = 0, k = start; j < len && w[j] == 0; j++, k++ ); + if ( j < len ) + return w[j] > 0 ? 1 : -1; + } + break; + case 2: + for ( j = 0, k = start; j < len && w[j] == 0; j++, k++ ); + if ( j < len ) + return w[j] > 0 ? 1 : -1; + break; + } + break; + } + } + return 0; +} + int cmpdl_matrix(int n,DL d1,DL d2) { int *v,*w,*t1,*t2; @@ -1525,8 +1621,8 @@ int cmpdl_matrix(int n,DL d1,DL d2) for ( i = 0, t1 = d1->d, t2 = d2->d, w = dp_dl_work; i < n; i++ ) w[i] = t1[i]-t2[i]; - len = dp_current_spec.ord.matrix.row; - matrix = dp_current_spec.ord.matrix.matrix; + len = dp_current_spec->ord.matrix.row; + matrix = dp_current_spec->ord.matrix.matrix; for ( j = 0; j < len; j++ ) { v = matrix[j]; for ( i = 0, s = 0; i < n; i++ ) @@ -1589,5 +1685,95 @@ DL remove_head_bucket(GeoBucket g,int nv) d = g->body[j]->body; g->body[j] = NEXT(g->body[j]); return d; + } +} + +/* DPV functions */ + +void adddv(VL vl,DPV p1,DPV p2,DPV *pr) +{ + int i,len; + DP *e; + + if ( !p1 || !p2 ) + error("adddv : invalid argument"); + else if ( p1->len != p2->len ) + error("adddv : size mismatch"); + else { + len = p1->len; + e = (DP *)MALLOC(p1->len*sizeof(DP)); + for ( i = 0; i < len; i++ ) + addd(vl,p1->body[i],p2->body[i],&e[i]); + MKDPV(len,e,*pr); + (*pr)->sugar = MAX(p1->sugar,p2->sugar); + } +} + +void subdv(VL vl,DPV p1,DPV p2,DPV *pr) +{ + int i,len; + DP *e; + + if ( !p1 || !p2 ) + error("subdv : invalid argument"); + else if ( p1->len != p2->len ) + error("subdv : size mismatch"); + else { + len = p1->len; + e = (DP *)MALLOC(p1->len*sizeof(DP)); + for ( i = 0; i < len; i++ ) + subd(vl,p1->body[i],p2->body[i],&e[i]); + MKDPV(len,e,*pr); + (*pr)->sugar = MAX(p1->sugar,p2->sugar); + } +} + +void chsgndv(DPV p1,DPV *pr) +{ + int i,len; + DP *e; + + if ( !p1 ) + error("subdv : invalid argument"); + else { + len = p1->len; + e = (DP *)MALLOC(p1->len*sizeof(DP)); + for ( i = 0; i < len; i++ ) + chsgnd(p1->body[i],&e[i]); + MKDPV(len,e,*pr); + (*pr)->sugar = p1->sugar; + } +} + +void muldv(VL vl,DP p1,DPV p2,DPV *pr) +{ + int i,len; + DP *e; + + len = p2->len; + e = (DP *)MALLOC(p2->len*sizeof(DP)); + if ( !p1 ) { + MKDPV(len,e,*pr); + (*pr)->sugar = 0; + } else { + for ( i = 0; i < len; i++ ) + muld(vl,p1,p2->body[i],&e[i]); + MKDPV(len,e,*pr); + (*pr)->sugar = p1->sugar + p2->sugar; + } +} + +int compdv(VL vl,DPV p1,DPV p2) +{ + int i,t,len; + + if ( p1->len != p2->len ) + error("compdv : size mismatch"); + else { + len = p1->len; + for ( i = 0; i < len; i++ ) + if ( t = compd(vl,p1->body[i],p2->body[i]) ) + return t; + return 0; } }