Annotation of OpenXM/src/asir-contrib/testing/noro/ndbf.rr, Revision 1.1
1.1 ! noro 1: /* requires 'primdec' */
! 2:
! 3: #define TMP_H hhhhhhhh
! 4: #define TMP_S ssssssss
! 5: #define TMP_DS dssssssss
! 6: #define TMP_T t
! 7: #define TMP_DT dt
! 8: #define TMP_Y1 yyyyyyyy1
! 9: #define TMP_DY1 dyyyyyyyy1
! 10: #define TMP_Y2 yyyyyyyy2
! 11: #define TMP_DY2 dyyyyyyyy2
! 12:
! 13: if (!module_definedp("gr")) load("gr")$ else{ }$
! 14: if (!module_definedp("primdec")) load("primdec")$ else{ }$
! 15: /* Empty for now. It will be used in a future. */
! 16:
! 17: /* toplevel */
! 18:
! 19: module ndbf$
! 20:
! 21: /* bfunction */
! 22:
! 23: localf bfunction, in_ww, in_ww_main, ann, ann_n$
! 24: localf ann0, psi, ww_weight, compare_first, generic_bfct$
! 25: localf generic_bfct_1, initial_part, bfct, indicial1, bfct_via_gbfct$
! 26: localf bfct_via_gbfct_weight, bfct_via_gbfct_weight_1, bfct_via_gbfct_weight_2$
! 27: localf weyl_minipolym, weyl_minipoly, weyl_nf, weyl_nf_quo_check$
! 28: localf weyl_nf_quo, weyl_nf_mod, b_subst, v_factorial, w_tdeg$
! 29: localf replace_vars_f, replace_vars_v, replace_var$
! 30: localf action_on_gfs, action_on_gfs_1$
! 31:
! 32: /* stratification */
! 33:
! 34: localf weyl_subst, bfactor, gen_a, gen_d$
! 35: localf gen_dm, elimination, weyl_ideal_quotient, psi0$
! 36: localf bf_strat, bf_strat_stage2, bf_strat_stage3, bf_local$
! 37: localf conv_tdt, merge_tower, stratify_bf, tdt_homogenize$
! 38: localf sing, tower_in_p, subst_vars, compute_exponent$
! 39: localf zero_inclusion, weyl_divide_by_right, elim_mat, int_cons$
! 40: localf ideal_intersection$
! 41:
! 42: def bfunction(F)
! 43: {
! 44: if ( type(Heu=getopt(heuristic)) == -1 ) Heu = 0;
! 45: if ( type(Vord=getopt(vord)) == -1 || type(Vord) != 4 ) Vord = 0;
! 46: if ( type(Wt=getopt(weight)) == -1 ) Wt = 0;
! 47: L = in_ww(F|weight=Wt,heuristic=Heu,vord=Vord);
! 48: Indata = L[0]; AllData = L[1]; VData = L[2];
! 49: GIN = Indata[0]; VDV = Indata[1]; WVDV = AllData[4];
! 50: W = Indata[4];
! 51: dp_set_weight(W);
! 52: B = weyl_minipoly(GIN,VDV,0,WVDV);
! 53: dp_set_weight(0);
! 54: return subst(B,s,-s-1);
! 55: }
! 56:
! 57: /*
! 58: returns [InData,AllData,VData]
! 59: InData = [GIN,VDV,V,DV,WtV]
! 60: AllData = [G0,GIN0,VDV0,W,WVDV,WtV0]
! 61: VData = [V0,DV0,T,DT]
! 62: GIN0 = ini_(-W,W)(G0)
! 63: WVDV = W[0]*V[0]*DV[0]+...
! 64: */
! 65:
! 66: def in_ww(F)
! 67: {
! 68: V = vars(F);
! 69: N = length(V);
! 70: D = newvect(N);
! 71: Wt = getopt(weight);
! 72: Vord = getopt(vord);
! 73: if ( type(Wt) != 4 ) {
! 74: if ( type(Vord) != 4 ) {
! 75: for ( I = 0; I < N; I++ )
! 76: D[I] = [deg(F,V[I]),V[I]];
! 77: qsort(D,compare_first);
! 78: for ( V = [], I = 0; I < N; I++ )
! 79: V = cons(D[I][1],V);
! 80: V = reverse(V);
! 81: } else
! 82: V = Vord;
! 83: for ( I = 0, Wt = []; I < N; I++ )
! 84: Wt = cons(1,Wt);
! 85: } else {
! 86: Wt1 = vector(N);
! 87: if ( type(Vord) != 4 ) {
! 88: for ( I = 0, F1 =F; I < N; I++ ) {
! 89: VI = Wt[2*I]; WI = Wt[2*I+1];
! 90: for ( J = 0; J < N; J++ )
! 91: if ( VI == V[J] ) break;
! 92: F1 = subst(F1,VI,VI^WI);
! 93: }
! 94: for ( I = 0; I < N; I++ )
! 95: D[I] = [deg(F1,V[I]),V[I]];
! 96: qsort(D,compare_first);
! 97: for ( V = [], I = 0; I < N; I++ )
! 98: V = cons(D[I][1],V);
! 99: V = reverse(V);
! 100: } else
! 101: V = Vord;
! 102: for ( I = 0; I < N; I++ ) {
! 103: VI = Wt[2*I]; WI = Wt[2*I+1];
! 104: for ( J = 0; J < N; J++ )
! 105: if ( VI == V[J] ) break;
! 106: Wt1[J] = WI;
! 107: }
! 108: Wt = vtol(Wt1);
! 109: }
! 110: Tdeg = w_tdeg(F,V,Wt);
! 111: /* weight for [t,x1,...,xn,dt,dx1,...,dxn,h] */
! 112: WtV = newvect(2*(N+1)+1);
! 113: WtV[0] = Tdeg; WtV[N+1] = 1; WtV[2*(N+1)] = 1;
! 114: /* wdeg(V[I])=Wt[I], wdeg(DV[I])=Tdeg-Wt[I]+1 */
! 115: for ( I = 1; I <= N; I++ ) {
! 116: WtV[I] = Wt[I-1]; WtV[N+1+I] = Tdeg-Wt[I-1]+1;
! 117: }
! 118: for ( I = N-1, DV = []; I >= 0; I-- )
! 119: DV = cons(strtov("d"+rtostr(V[I])),DV);
! 120:
! 121: B = [TMP_T-F];
! 122: for ( I = 0; I < N; I++ ) {
! 123: B = cons(DV[I]+diff(F,V[I])*TMP_DT,B);
! 124: }
! 125: V1 = cons(TMP_T,V); DV1 = cons(TMP_DT,DV);
! 126: W = newvect(N+1); W[0] = 1;
! 127: VW1 = [V1,DV1,WtV,W];
! 128:
! 129: /* weight for [x1,...,xn,t,dx1,...,dxn,dt,h] */
! 130: WtV = newvect(2*(N+1)+1);
! 131: WtV[N] = Tdeg; WtV[2*N+1] = 1; WtV[2*(N+1)] = 1;
! 132: for ( I = 0; I < N; I++ ) {
! 133: WtV[I] = Wt[I]; WtV[N+1+I] = Tdeg-Wt[I]+1;
! 134: }
! 135: for ( I = N-1, DV = []; I >= 0; I-- )
! 136: DV = cons(strtov("d"+rtostr(V[I])),DV);
! 137:
! 138: B = [TMP_T-F];
! 139: for ( I = 0; I < N; I++ ) {
! 140: B = cons(DV[I]+diff(F,V[I])*TMP_DT,B);
! 141: }
! 142: V1 = append(V,[TMP_T]); DV1 = append(DV,[TMP_DT]);
! 143: W = newvect(N+1); W[N] = 1;
! 144: VW2 = [V1,DV1,WtV,W];
! 145:
! 146: Heu = getopt(heuristic);
! 147: if ( type(Heu) != -1 && Heu )
! 148: L = in_ww_main(B,VW1,VW2);
! 149: else
! 150: L = in_ww_main(B,VW1,0);
! 151: return append(L,[[V,DV,TMP_T,TMP_DT,Wt]]);
! 152: }
! 153:
! 154: /*
! 155: returns [InData,AllData]
! 156: InData = [GIN,VDV,V,DV,WtV]
! 157: AllData = [G0,GIN0,VDV0,W,WVDV,WtV0]
! 158: GIN0 = ini_(-W,W)(G0)
! 159: WVDV = W[0]*V[0]*DV[0]+...
! 160: */
! 161:
! 162: def in_ww_main(F,VW1,VW2)
! 163: {
! 164: V = VW1[0]; DV = VW1[1]; WtV = VW1[2]; W = VW1[3];
! 165: dp_set_weight(WtV);
! 166:
! 167: N = length(V);
! 168: N2 = N*2;
! 169:
! 170: /* If W is a list, convert it to a vector */
! 171: if ( type(W) == 4 )
! 172: W = newvect(length(W),W);
! 173: dp_weyl_set_weight(W);
! 174:
! 175: /* create a term order M in D<x,d> (DRL) */
! 176: M = newmat(N2,N2);
! 177: for ( J = 0; J < N2; J++ )
! 178: M[0][J] = 1;
! 179: for ( I = 1; I < N2; I++ )
! 180: M[I][N2-I] = -1;
! 181:
! 182: VDV = append(V,DV);
! 183:
! 184: /* create a non-term order MW in D<x,d> */
! 185: MW = newmat(N2+1,N2);
! 186: for ( J = 0; J < N; J++ )
! 187: MW[0][J] = -W[J];
! 188: for ( ; J < N2; J++ )
! 189: MW[0][J] = W[J-N];
! 190: for ( I = 1; I <= N2; I++ )
! 191: for ( J = 0; J < N2; J++ )
! 192: MW[I][J] = M[I-1][J];
! 193:
! 194: /* create a homogenized term order MWH in D<x,d,h> */
! 195: MWH = newmat(N2+2,N2+1);
! 196: for ( J = 0; J <= N2; J++ )
! 197: MWH[0][J] = 1;
! 198: for ( I = 1; I <= N2+1; I++ )
! 199: for ( J = 0; J < N2; J++ )
! 200: MWH[I][J] = MW[I-1][J];
! 201:
! 202: /* homogenize F */
! 203: VDVH = append(VDV,[TMP_H]);
! 204: FH = map(dp_dtop,map(dp_homo,map(dp_ptod,F,VDV)),VDVH);
! 205:
! 206: /* compute a groebner basis of FH w.r.t. MWH */
! 207: /* dp_gr_flags(["Top",1,"NoRA",1]); */
! 208: #if 0
! 209: GH = dp_weyl_gr_main(FH,VDVH,0,0,11);
! 210: #else
! 211: #if 0
! 212: GH = dp_weyl_gr_main(FH,VDVH,0,-1,11);
! 213: #else
! 214: GH = nd_weyl_gr_trace(FH,VDVH,0,-1,11);
! 215: #endif
! 216: #endif
! 217: /* dp_gr_flags(["Top",0,"NoRA",0]); */
! 218:
! 219: /* dehomigenize GH */
! 220: G = map(subst,GH,TMP_H,1);
! 221:
! 222: /* G is a groebner basis w.r.t. a non term order MW */
! 223: /* take the initial part w.r.t. (-W,W) */
! 224: GIN = map(initial_part,G,VDV,MW,W);
! 225:
! 226: /* GIN is a groebner basis w.r.t. a term order M */
! 227: /* As -W+W=0, gr_(-W,W)(D<x,d>) = D<x,d> */
! 228:
! 229: /* find b(W1*x1*d1+...+WN*xN*dN) in Id(GIN) */
! 230: for ( I = 0, T = 0; I < N; I++ )
! 231: T += W[I]*V[I]*DV[I];
! 232:
! 233: AllData = [G,GIN,VDV,W,T,WtV];
! 234: if ( VW2 ) {
! 235: V2 = VW2[0]; DV2 = VW2[1]; WtV2 = VW2[2];
! 236: VDV2 = append(V2,DV2);
! 237: dp_set_weight(WtV2);
! 238: GIN2 = nd_weyl_gr_trace(GIN,VDV2,0,-1,0);
! 239: InData = [GIN2,VDV2,V2,DV2,WtV2];
! 240: } else {
! 241: if ( 0 ) {
! 242: dp_set_weight(WtV);
! 243: GIN1 = nd_weyl_gr_postproc(GIN,VDV,0,0,0);
! 244: InData = [GIN1,VDV,V,DV,WtV];
! 245: } else
! 246: InData = [GIN,VDV,V,DV,WtV];
! 247: }
! 248:
! 249: /* B = weyl_minipoly(GIN2,VDV2,0,T); */ /* M represents DRL order */
! 250: WtV = dp_set_weight();
! 251: dp_set_weight(0);
! 252:
! 253: return [InData,AllData];
! 254: }
! 255:
! 256: /* annihilating ideal of F^s */
! 257:
! 258: def ann(F)
! 259: {
! 260: if ( member(s,vars(F)) )
! 261: error("ann : the variable 's' is reserved.");
! 262: V = vars(F);
! 263: N = length(V);
! 264: D = newvect(N);
! 265:
! 266: for ( I = 0; I < N; I++ )
! 267: D[I] = [deg(F,V[I]),V[I]];
! 268: qsort(D,compare_first);
! 269: for ( V = [], I = N-1; I >= 0; I-- )
! 270: V = cons(D[I][1],V);
! 271:
! 272: for ( I = N-1, DV = []; I >= 0; I-- )
! 273: DV = cons(strtov("d"+rtostr(V[I])),DV);
! 274:
! 275: W = append([TMP_Y1,TMP_Y2,TMP_T],V);
! 276: DW = append([TMP_DY1,TMP_DY2,TMP_DT],DV);
! 277:
! 278: B = [1-TMP_Y1*TMP_Y2,TMP_T-TMP_Y1*F];
! 279: for ( I = 0; I < N; I++ ) {
! 280: B = cons(DV[I]+TMP_Y1*diff(F,V[I])*TMP_DT,B);
! 281: }
! 282:
! 283: /* homogenized (heuristics) */
! 284: dp_nelim(2);
! 285: G0 = nd_weyl_gr(B,append(W,DW),0,[[0,2],[0,length(W)*2-2]]);
! 286: G1 = [];
! 287: for ( T = G0; T != []; T = cdr(T) ) {
! 288: E = car(T); VL = vars(E);
! 289: if ( !member(TMP_Y1,VL) && !member(TMP_Y2,VL) )
! 290: G1 = cons(E,G1);
! 291: }
! 292: G2 = map(psi,G1,TMP_T,TMP_DT);
! 293: G3 = map(subst,G2,TMP_T,-1-s);
! 294: return G3;
! 295: }
! 296:
! 297: /* F = [F0,F1,...] */
! 298:
! 299: def ann_n(F)
! 300: {
! 301: L = length(F);
! 302: V = vars(F);
! 303: N = length(V);
! 304: D = newvect(N);
! 305:
! 306: for ( I = N-1, DV = []; I >= 0; I-- )
! 307: DV = cons(strtov("d"+rtostr(V[I])),DV);
! 308: W = []; DW = [];
! 309: for ( I = L-1; I >= 0; I-- ) {
! 310: SI = rtostr(I);
! 311: W = cons(strtov("t"+SI),W);
! 312: DW = cons(strtov("dt"+SI),DW);
! 313: }
! 314: U = []; DU = [];
! 315: for ( I = L-1; I >= 0; I-- ) {
! 316: SI = rtostr(I);
! 317: U = append([strtov("u"+SI),strtov("v"+SI)],U);
! 318: DU = append([strtov("du"+SI),strtov("dv"+SI)],DU);
! 319: }
! 320:
! 321: B = [];
! 322: for ( I = 0; I < N; I++ ) {
! 323: T = DV[I];
! 324: for ( J = 0; J < L; J++ )
! 325: T += U[2*J]*diff(F[J],V[I])*DW[J];
! 326: B = cons(T,B);
! 327: }
! 328: for ( I = 0; I < L; I++ )
! 329: B = append([W[I]-U[2*I]*F[I],1-U[2*I]*U[2*I+1]],B);
! 330:
! 331: VA = append(U,append(W,V));
! 332: DVA = append(DU,append(DW,DV));
! 333: VDV = append(VA,DVA);
! 334: G0 = nd_weyl_gr(B,VDV,0,[[0,2*L],[0,length(VDV)-2*L]]);
! 335: G1 = [];
! 336: for ( T = G0; T != []; T = cdr(T) ) {
! 337: E = car(T); VL = vars(E);
! 338: for ( TV = U; TV != []; TV = cdr(TV) )
! 339: if ( member(car(TV),VL) ) break;
! 340: if ( TV == [] )
! 341: G1 = cons(E,G1);
! 342: }
! 343: G2 = G1;
! 344: for ( I = 0; I < L; I++ ) {
! 345: G2 = map(psi,G2,W[I],DW[I]);
! 346: G2 = map(subst,G2,W[I],-1-strtov("s"+rtostr(I)));
! 347: }
! 348: return G2;
! 349: }
! 350:
! 351: /*
! 352: * compute J_f|s=r, where r = the minimal integral root of global b_f(s)
! 353: * ann0(F) returns [MinRoot,Ideal]
! 354: */
! 355:
! 356: def ann0(F)
! 357: {
! 358: F = subst(F,s,TMP_S);
! 359: Ann = ann(F);
! 360: Bf = bfunction(F);
! 361:
! 362: FList = cdr(fctr(Bf));
! 363: for ( T = FList, Min = 0; T != []; T = cdr(T) ) {
! 364: LF = car(car(T));
! 365: Root = -coef(LF,0)/coef(LF,1);
! 366: if ( dn(Root) == 1 && Root < Min )
! 367: Min = Root;
! 368: }
! 369: return [Min,map(subst,Ann,s,Min,TMP_S,s,TMP_DS,ds)];
! 370: }
! 371:
! 372: def psi0(F,T,DT)
! 373: {
! 374: D = dp_ptod(F,[T,DT]);
! 375: Wmax = ww_weight(D);
! 376: D1 = dp_rest(D);
! 377: for ( ; D1; D1 = dp_rest(D1) )
! 378: if ( ww_weight(D1) > Wmax )
! 379: Wmax = ww_weight(D1);
! 380: for ( D1 = D, Dmax = 0; D1; D1 = dp_rest(D1) )
! 381: if ( ww_weight(D1) == Wmax )
! 382: Dmax += dp_hm(D1);
! 383: if ( Wmax >= 0 )
! 384: Dmax = dp_weyl_mul(<<Wmax,0>>,Dmax);
! 385: else
! 386: Dmax = dp_weyl_mul(<<0,-Wmax>>,Dmax);
! 387: Rmax = dp_dtop(Dmax,[T,DT]);
! 388: return Rmax;
! 389: }
! 390:
! 391: def psi(F,T,DT)
! 392: {
! 393: Rmax = psi0(F,T,DT);
! 394: R = b_subst(subst(Rmax,DT,1),T);
! 395: return R;
! 396: }
! 397:
! 398: def ww_weight(D)
! 399: {
! 400: V = dp_etov(D);
! 401: return V[1]-V[0];
! 402: }
! 403:
! 404: def compare_first(A,B)
! 405: {
! 406: A0 = car(A);
! 407: B0 = car(B);
! 408: if ( A0 > B0 )
! 409: return 1;
! 410: else if ( A0 < B0 )
! 411: return -1;
! 412: else
! 413: return 0;
! 414: }
! 415:
! 416: /* generic b-function w.r.t. weight vector W */
! 417:
! 418: def generic_bfct(F,V,DV,W)
! 419: {
! 420: N = length(V);
! 421: N2 = N*2;
! 422:
! 423: /* If W is a list, convert it to a vector */
! 424: if ( type(W) == 4 )
! 425: W = newvect(length(W),W);
! 426: dp_weyl_set_weight(W);
! 427:
! 428: /* create a term order M in D<x,d> (DRL) */
! 429: M = newmat(N2,N2);
! 430: for ( J = 0; J < N2; J++ )
! 431: M[0][J] = 1;
! 432: for ( I = 1; I < N2; I++ )
! 433: M[I][N2-I] = -1;
! 434:
! 435: VDV = append(V,DV);
! 436:
! 437: /* create a non-term order MW in D<x,d> */
! 438: MW = newmat(N2+1,N2);
! 439: for ( J = 0; J < N; J++ )
! 440: MW[0][J] = -W[J];
! 441: for ( ; J < N2; J++ )
! 442: MW[0][J] = W[J-N];
! 443: for ( I = 1; I <= N2; I++ )
! 444: for ( J = 0; J < N2; J++ )
! 445: MW[I][J] = M[I-1][J];
! 446:
! 447: /* create a homogenized term order MWH in D<x,d,h> */
! 448: MWH = newmat(N2+2,N2+1);
! 449: for ( J = 0; J <= N2; J++ )
! 450: MWH[0][J] = 1;
! 451: for ( I = 1; I <= N2+1; I++ )
! 452: for ( J = 0; J < N2; J++ )
! 453: MWH[I][J] = MW[I-1][J];
! 454:
! 455: /* homogenize F */
! 456: VDVH = append(VDV,[TMP_H]);
! 457: FH = map(dp_dtop,map(dp_homo,map(dp_ptod,F,VDV)),VDVH);
! 458:
! 459: /* compute a groebner basis of FH w.r.t. MWH */
! 460: dp_gr_flags(["Top",1,"NoRA",1]);
! 461: GH = dp_weyl_gr_main(FH,VDVH,0,1,11);
! 462: dp_gr_flags(["Top",0,"NoRA",0]);
! 463:
! 464: /* dehomigenize GH */
! 465: G = map(subst,GH,TMP_H,1);
! 466:
! 467: /* G is a groebner basis w.r.t. a non term order MW */
! 468: /* take the initial part w.r.t. (-W,W) */
! 469: GIN = map(initial_part,G,VDV,MW,W);
! 470:
! 471: /* GIN is a groebner basis w.r.t. a term order M */
! 472: /* As -W+W=0, gr_(-W,W)(D<x,d>) = D<x,d> */
! 473:
! 474: /* find b(W1*x1*d1+...+WN*xN*dN) in Id(GIN) */
! 475: for ( I = 0, T = 0; I < N; I++ )
! 476: T += W[I]*V[I]*DV[I];
! 477: B = weyl_minipoly(GIN,VDV,0,T); /* M represents DRL order */
! 478: return B;
! 479: }
! 480:
! 481: /* all term reduction + interreduce */
! 482: def generic_bfct_1(F,V,DV,W)
! 483: {
! 484: N = length(V);
! 485: N2 = N*2;
! 486:
! 487: /* If W is a list, convert it to a vector */
! 488: if ( type(W) == 4 )
! 489: W = newvect(length(W),W);
! 490: dp_weyl_set_weight(W);
! 491:
! 492: /* create a term order M in D<x,d> (DRL) */
! 493: M = newmat(N2,N2);
! 494: for ( J = 0; J < N2; J++ )
! 495: M[0][J] = 1;
! 496: for ( I = 1; I < N2; I++ )
! 497: M[I][N2-I] = -1;
! 498:
! 499: VDV = append(V,DV);
! 500:
! 501: /* create a non-term order MW in D<x,d> */
! 502: MW = newmat(N2+1,N2);
! 503: for ( J = 0; J < N; J++ )
! 504: MW[0][J] = -W[J];
! 505: for ( ; J < N2; J++ )
! 506: MW[0][J] = W[J-N];
! 507: for ( I = 1; I <= N2; I++ )
! 508: for ( J = 0; J < N2; J++ )
! 509: MW[I][J] = M[I-1][J];
! 510:
! 511: /* create a homogenized term order MWH in D<x,d,h> */
! 512: MWH = newmat(N2+2,N2+1);
! 513: for ( J = 0; J <= N2; J++ )
! 514: MWH[0][J] = 1;
! 515: for ( I = 1; I <= N2+1; I++ )
! 516: for ( J = 0; J < N2; J++ )
! 517: MWH[I][J] = MW[I-1][J];
! 518:
! 519: /* homogenize F */
! 520: VDVH = append(VDV,[TMP_H]);
! 521: FH = map(dp_dtop,map(dp_homo,map(dp_ptod,F,VDV)),VDVH);
! 522:
! 523: /* compute a groebner basis of FH w.r.t. MWH */
! 524: /* dp_gr_flags(["Top",1,"NoRA",1]); */
! 525: GH = dp_weyl_gr_main(FH,VDVH,0,1,11);
! 526: /* dp_gr_flags(["Top",0,"NoRA",0]); */
! 527:
! 528: /* dehomigenize GH */
! 529: G = map(subst,GH,TMP_H,1);
! 530:
! 531: /* G is a groebner basis w.r.t. a non term order MW */
! 532: /* take the initial part w.r.t. (-W,W) */
! 533: GIN = map(initial_part,G,VDV,MW,W);
! 534:
! 535: /* GIN is a groebner basis w.r.t. a term order M */
! 536: /* As -W+W=0, gr_(-W,W)(D<x,d>) = D<x,d> */
! 537:
! 538: /* find b(W1*x1*d1+...+WN*xN*dN) in Id(GIN) */
! 539: for ( I = 0, T = 0; I < N; I++ )
! 540: T += W[I]*V[I]*DV[I];
! 541: B = weyl_minipoly(GIN,VDV,0,T); /* M represents DRL order */
! 542: return B;
! 543: }
! 544:
! 545: def initial_part(F,V,MW,W)
! 546: {
! 547: N2 = length(V);
! 548: N = N2/2;
! 549: dp_ord(MW);
! 550: DF = dp_ptod(F,V);
! 551: R = dp_hm(DF);
! 552: DF = dp_rest(DF);
! 553:
! 554: E = dp_etov(R);
! 555: for ( I = 0, TW = 0; I < N; I++ )
! 556: TW += W[I]*(-E[I]+E[N+I]);
! 557: RW = TW;
! 558:
! 559: for ( ; DF; DF = dp_rest(DF) ) {
! 560: E = dp_etov(DF);
! 561: for ( I = 0, TW = 0; I < N; I++ )
! 562: TW += W[I]*(-E[I]+E[N+I]);
! 563: if ( TW == RW )
! 564: R += dp_hm(DF);
! 565: else if ( TW < RW )
! 566: break;
! 567: else
! 568: error("initial_part : cannot happen");
! 569: }
! 570: return dp_dtop(R,V);
! 571:
! 572: }
! 573:
! 574: /* b-function of F ? */
! 575:
! 576: def bfct(F)
! 577: {
! 578: /* XXX */
! 579: F = replace_vars_f(F);
! 580:
! 581: V = vars(F);
! 582: N = length(V);
! 583: D = newvect(N);
! 584:
! 585: for ( I = 0; I < N; I++ )
! 586: D[I] = [deg(F,V[I]),V[I]];
! 587: qsort(D,compare_first);
! 588: for ( V = [], I = 0; I < N; I++ )
! 589: V = cons(D[I][1],V);
! 590: for ( I = N-1, DV = []; I >= 0; I-- )
! 591: DV = cons(strtov("d"+rtostr(V[I])),DV);
! 592: V1 = cons(s,V); DV1 = cons(ds,DV);
! 593:
! 594: G0 = indicial1(F,reverse(V));
! 595: G1 = dp_weyl_gr_main(G0,append(V1,DV1),0,1,0);
! 596: Minipoly = weyl_minipoly(G1,append(V1,DV1),0,s);
! 597: return Minipoly;
! 598: }
! 599:
! 600: /* called from bfct() only */
! 601:
! 602: def indicial1(F,V)
! 603: {
! 604: W = append([y1,t],V);
! 605: N = length(V);
! 606: B = [t-y1*F];
! 607: for ( I = N-1, DV = []; I >= 0; I-- )
! 608: DV = cons(strtov("d"+rtostr(V[I])),DV);
! 609: DW = append([dy1,dt],DV);
! 610: for ( I = 0; I < N; I++ ) {
! 611: B = cons(DV[I]+y1*diff(F,V[I])*dt,B);
! 612: }
! 613: dp_nelim(1);
! 614:
! 615: /* homogenized (heuristics) */
! 616: G0 = dp_weyl_gr_main(B,append(W,DW),1,0,6);
! 617: G1 = map(subst,G0,y1,1);
! 618: G2 = map(psi,G1,t,dt);
! 619: G3 = map(subst,G2,t,-s-1);
! 620: return G3;
! 621: }
! 622:
! 623: /* b-function computation via generic_bfct() (experimental) */
! 624:
! 625: def bfct_via_gbfct(F)
! 626: {
! 627: V = vars(F);
! 628: N = length(V);
! 629: D = newvect(N);
! 630:
! 631: for ( I = 0; I < N; I++ )
! 632: D[I] = [deg(F,V[I]),V[I]];
! 633: qsort(D,compare_first);
! 634: for ( V = [], I = 0; I < N; I++ )
! 635: V = cons(D[I][1],V);
! 636: V = reverse(V);
! 637: for ( I = N-1, DV = []; I >= 0; I-- )
! 638: DV = cons(strtov("d"+rtostr(V[I])),DV);
! 639:
! 640: B = [TMP_T-F];
! 641: for ( I = 0; I < N; I++ ) {
! 642: B = cons(DV[I]+diff(F,V[I])*TMP_DT,B);
! 643: }
! 644: V1 = cons(TMP_T,V); DV1 = cons(TMP_DT,DV);
! 645: W = newvect(N+1);
! 646: W[0] = 1;
! 647: R = generic_bfct(B,V1,DV1,W);
! 648:
! 649: return subst(R,s,-s-1);
! 650: }
! 651:
! 652: /* use an order s.t. [t,x,y,z,...,dt,dx,dy,dz,...,h] */
! 653:
! 654: def bfct_via_gbfct_weight(F,V)
! 655: {
! 656: N = length(V);
! 657: D = newvect(N);
! 658: Wt = getopt(weight);
! 659: if ( type(Wt) != 4 ) {
! 660: for ( I = 0, Wt = []; I < N; I++ )
! 661: Wt = cons(1,Wt);
! 662: }
! 663: Tdeg = w_tdeg(F,V,Wt);
! 664: WtV = newvect(2*(N+1)+1);
! 665: WtV[0] = Tdeg;
! 666: WtV[N+1] = 1;
! 667: /* wdeg(V[I])=Wt[I], wdeg(DV[I])=Tdeg-Wt[I]+1 */
! 668: for ( I = 1; I <= N; I++ ) {
! 669: WtV[I] = Wt[I-1];
! 670: WtV[N+1+I] = Tdeg-Wt[I-1]+1;
! 671: }
! 672: WtV[2*(N+1)] = 1;
! 673: dp_set_weight(WtV);
! 674: for ( I = N-1, DV = []; I >= 0; I-- )
! 675: DV = cons(strtov("d"+rtostr(V[I])),DV);
! 676:
! 677: B = [TMP_T-F];
! 678: for ( I = 0; I < N; I++ ) {
! 679: B = cons(DV[I]+diff(F,V[I])*TMP_DT,B);
! 680: }
! 681: V1 = cons(TMP_T,V); DV1 = cons(TMP_DT,DV);
! 682: W = newvect(N+1);
! 683: W[0] = 1;
! 684: R = generic_bfct_1(B,V1,DV1,W);
! 685: dp_set_weight(0);
! 686: return subst(R,s,-s-1);
! 687: }
! 688:
! 689: /* use an order s.t. [x,y,z,...,t,dx,dy,dz,...,dt,h] */
! 690:
! 691: def bfct_via_gbfct_weight_1(F,V)
! 692: {
! 693: N = length(V);
! 694: D = newvect(N);
! 695: Wt = getopt(weight);
! 696: if ( type(Wt) != 4 ) {
! 697: for ( I = 0, Wt = []; I < N; I++ )
! 698: Wt = cons(1,Wt);
! 699: }
! 700: Tdeg = w_tdeg(F,V,Wt);
! 701: WtV = newvect(2*(N+1)+1);
! 702: /* wdeg(V[I])=Wt[I], wdeg(DV[I])=Tdeg-Wt[I]+1 */
! 703: for ( I = 0; I < N; I++ ) {
! 704: WtV[I] = Wt[I];
! 705: WtV[N+1+I] = Tdeg-Wt[I]+1;
! 706: }
! 707: WtV[N] = Tdeg;
! 708: WtV[2*N+1] = 1;
! 709: WtV[2*(N+1)] = 1;
! 710: dp_set_weight(WtV);
! 711: for ( I = N-1, DV = []; I >= 0; I-- )
! 712: DV = cons(strtov("d"+rtostr(V[I])),DV);
! 713:
! 714: B = [TMP_T-F];
! 715: for ( I = 0; I < N; I++ ) {
! 716: B = cons(DV[I]+diff(F,V[I])*TMP_DT,B);
! 717: }
! 718: V1 = append(V,[TMP_T]); DV1 = append(DV,[TMP_DT]);
! 719: W = newvect(N+1);
! 720: W[N] = 1;
! 721: R = generic_bfct_1(B,V1,DV1,W);
! 722: dp_set_weight(0);
! 723: return subst(R,s,-s-1);
! 724: }
! 725:
! 726: def bfct_via_gbfct_weight_2(F,V)
! 727: {
! 728: N = length(V);
! 729: D = newvect(N);
! 730: Wt = getopt(weight);
! 731: if ( type(Wt) != 4 ) {
! 732: for ( I = 0, Wt = []; I < N; I++ )
! 733: Wt = cons(1,Wt);
! 734: }
! 735: Tdeg = w_tdeg(F,V,Wt);
! 736:
! 737: /* a weight for the first GB computation */
! 738: /* [t,x1,...,xn,dt,dx1,...,dxn,h] */
! 739: WtV = newvect(2*(N+1)+1);
! 740: WtV[0] = Tdeg;
! 741: WtV[N+1] = 1;
! 742: WtV[2*(N+1)] = 1;
! 743: /* wdeg(V[I])=Wt[I], wdeg(DV[I])=Tdeg-Wt[I]+1 */
! 744: for ( I = 1; I <= N; I++ ) {
! 745: WtV[I] = Wt[I-1];
! 746: WtV[N+1+I] = Tdeg-Wt[I-1]+1;
! 747: }
! 748: dp_set_weight(WtV);
! 749:
! 750: /* a weight for the second GB computation */
! 751: /* [x1,...,xn,t,dx1,...,dxn,dt,h] */
! 752: WtV2 = newvect(2*(N+1)+1);
! 753: WtV2[N] = Tdeg;
! 754: WtV2[2*N+1] = 1;
! 755: WtV2[2*(N+1)] = 1;
! 756: for ( I = 0; I < N; I++ ) {
! 757: WtV2[I] = Wt[I];
! 758: WtV2[N+1+I] = Tdeg-Wt[I]+1;
! 759: }
! 760:
! 761: for ( I = N-1, DV = []; I >= 0; I-- )
! 762: DV = cons(strtov("d"+rtostr(V[I])),DV);
! 763:
! 764: B = [TMP_T-F];
! 765: for ( I = 0; I < N; I++ ) {
! 766: B = cons(DV[I]+diff(F,V[I])*TMP_DT,B);
! 767: }
! 768: V1 = cons(TMP_T,V); DV1 = cons(TMP_DT,DV);
! 769: V2 = append(V,[TMP_T]); DV2 = append(DV,[TMP_DT]);
! 770: W = newvect(N+1,[1]);
! 771: dp_weyl_set_weight(W);
! 772:
! 773: VDV = append(V1,DV1);
! 774: N1 = length(V1);
! 775: N2 = N1*2;
! 776:
! 777: /* create a non-term order MW in D<x,d> */
! 778: MW = newmat(N2+1,N2);
! 779: for ( J = 0; J < N1; J++ ) {
! 780: MW[0][J] = -W[J]; MW[0][N1+J] = W[J];
! 781: }
! 782: for ( J = 0; J < N2; J++ ) MW[1][J] = 1;
! 783: for ( I = 2; I <= N2; I++ ) MW[I][N2-I+1] = -1;
! 784:
! 785: /* homogenize F */
! 786: VDVH = append(VDV,[TMP_H]);
! 787: FH = map(dp_dtop,map(dp_homo,map(dp_ptod,B,VDV)),VDVH);
! 788:
! 789: /* compute a groebner basis of FH w.r.t. MWH */
! 790: GH = dp_weyl_gr_main(FH,VDVH,0,1,11);
! 791:
! 792: /* dehomigenize GH */
! 793: G = map(subst,GH,TMP_H,1);
! 794:
! 795: /* G is a groebner basis w.r.t. a non term order MW */
! 796: /* take the initial part w.r.t. (-W,W) */
! 797: GIN = map(initial_part,G,VDV,MW,W);
! 798:
! 799: /* GIN is a groebner basis w.r.t. a term order M */
! 800: /* As -W+W=0, gr_(-W,W)(D<x,d>) = D<x,d> */
! 801:
! 802: /* find b(W1*x1*d1+...+WN*xN*dN) in Id(GIN) */
! 803: for ( I = 0, T = 0; I < N1; I++ )
! 804: T += W[I]*V1[I]*DV1[I];
! 805:
! 806: /* change of ordering from VDV to VDV2 */
! 807: VDV2 = append(V2,DV2);
! 808: dp_set_weight(WtV2);
! 809: for ( Pind = 0; ; Pind++ ) {
! 810: Prime = lprime(Pind);
! 811: GIN2 = dp_weyl_gr_main(GIN,VDV2,0,-Prime,0);
! 812: if ( GIN2 ) break;
! 813: }
! 814:
! 815: R = weyl_minipoly(GIN2,VDV2,0,T); /* M represents DRL order */
! 816: dp_set_weight(0);
! 817: return subst(R,s,-s-1);
! 818: }
! 819:
! 820: /* minimal polynomial of s; modular computation */
! 821:
! 822: def weyl_minipolym(G,V,O,M,V0)
! 823: {
! 824: N = length(V);
! 825: Len = length(G);
! 826: dp_ord(O);
! 827: setmod(M);
! 828: PS = newvect(Len);
! 829: PS0 = newvect(Len);
! 830:
! 831: for ( I = 0, T = G; T != []; T = cdr(T), I++ )
! 832: PS0[I] = dp_ptod(car(T),V);
! 833: for ( I = 0, T = G; T != []; T = cdr(T), I++ )
! 834: PS[I] = dp_mod(dp_ptod(car(T),V),M,[]);
! 835:
! 836: for ( I = Len - 1, GI = []; I >= 0; I-- )
! 837: GI = cons(I,GI);
! 838:
! 839: U = dp_mod(dp_ptod(V0,V),M,[]);
! 840: U = dp_weyl_nf_mod(GI,U,PS,1,M);
! 841:
! 842: T = dp_mod(<<0>>,M,[]);
! 843: TT = dp_mod(dp_ptod(1,V),M,[]);
! 844: G = H = [[TT,T]];
! 845:
! 846: for ( I = 1; ; I++ ) {
! 847: if ( dp_gr_print() )
! 848: print(".",2);
! 849: T = dp_mod(<<I>>,M,[]);
! 850:
! 851: TT = dp_weyl_nf_mod(GI,dp_weyl_mul_mod(TT,U,M),PS,1,M);
! 852: H = cons([TT,T],H);
! 853: L = dp_lnf_mod([TT,T],G,M);
! 854: if ( !L[0] ) {
! 855: if ( dp_gr_print() )
! 856: print("");
! 857: return dp_dtop(L[1],[t]); /* XXX */
! 858: } else
! 859: G = insert(G,L);
! 860: }
! 861: }
! 862:
! 863: /* minimal polynomial of s over Q */
! 864:
! 865: def weyl_minipoly(G0,V0,O0,P)
! 866: {
! 867: HM = hmlist(G0,V0,O0);
! 868:
! 869: N = length(V0);
! 870: Len = length(G0);
! 871: dp_ord(O0);
! 872: PS = newvect(Len);
! 873: for ( I = 0, T = G0, HL = []; T != []; T = cdr(T), I++ )
! 874: PS[I] = dp_ptod(car(T),V0);
! 875: for ( I = Len - 1, GI = []; I >= 0; I-- )
! 876: GI = cons(I,GI);
! 877: PSM = newvect(Len);
! 878: DP = dp_ptod(P,V0);
! 879:
! 880: for ( Pind = 0; ; Pind++ ) {
! 881: Prime = lprime(Pind);
! 882: if ( !valid_modulus(HM,Prime) )
! 883: continue;
! 884: setmod(Prime);
! 885: for ( I = 0, T = G0, HL = []; T != []; T = cdr(T), I++ )
! 886: PSM[I] = dp_mod(dp_ptod(car(T),V0),Prime,[]);
! 887:
! 888: NFP = weyl_nf(GI,DP,1,PS);
! 889: NFPM = dp_mod(NFP[0],Prime,[])/ptomp(NFP[1],Prime);
! 890:
! 891: NF = [[dp_ptod(1,V0),1]];
! 892: LCM = 1;
! 893:
! 894: TM = dp_mod(<<0>>,Prime,[]);
! 895: TTM = dp_mod(dp_ptod(1,V0),Prime,[]);
! 896: GM = NFM = [[TTM,TM]];
! 897:
! 898: for ( D = 1; ; D++ ) {
! 899: if ( dp_gr_print() )
! 900: print(".",2);
! 901: NFPrev = car(NF);
! 902: NFJ = weyl_nf(GI,
! 903: dp_weyl_mul(NFP[0],NFPrev[0]),NFP[1]*NFPrev[1],PS);
! 904: NFJ = remove_cont(NFJ);
! 905: NF = cons(NFJ,NF);
! 906: LCM = ilcm(LCM,NFJ[1]);
! 907:
! 908: /* modular computation */
! 909: TM = dp_mod(<<D>>,Prime,[]);
! 910: TTM = dp_mod(NFJ[0],Prime,[])/ptomp(NFJ[1],Prime);
! 911: NFM = cons([TTM,TM],NFM);
! 912: LM = dp_lnf_mod([TTM,TM],GM,Prime);
! 913: if ( !LM[0] )
! 914: break;
! 915: else
! 916: GM = insert(GM,LM);
! 917: }
! 918:
! 919: if ( dp_gr_print() )
! 920: print("");
! 921: U = NF[0][0]*idiv(LCM,NF[0][1]);
! 922: Coef = [];
! 923: for ( J = D-1; J >= 0; J-- ) {
! 924: Coef = cons(strtov("u"+rtostr(J)),Coef);
! 925: U += car(Coef)*NF[D-J][0]*idiv(LCM,NF[D-J][1]);
! 926: }
! 927:
! 928: for ( UU = U, Eq = []; UU; UU = dp_rest(UU) )
! 929: Eq = cons(dp_hc(UU),Eq);
! 930: M = etom([Eq,Coef]);
! 931: B = henleq(M,Prime);
! 932: if ( dp_gr_print() )
! 933: print("");
! 934: if ( B ) {
! 935: R = 0;
! 936: for ( I = 0; I < D; I++ )
! 937: R += B[0][I]*s^I;
! 938: R += B[1]*s^D;
! 939: return R;
! 940: }
! 941: }
! 942: }
! 943:
! 944: def weyl_nf(B,G,M,PS)
! 945: {
! 946: for ( D = 0; G; ) {
! 947: for ( U = 0, L = B; L != []; L = cdr(L) ) {
! 948: if ( dp_redble(G,R=PS[car(L)]) > 0 ) {
! 949: GCD = igcd(dp_hc(G),dp_hc(R));
! 950: CG = idiv(dp_hc(R),GCD); CR = idiv(dp_hc(G),GCD);
! 951: U = CG*G-dp_weyl_mul(CR*dp_subd(G,R),R);
! 952: if ( !U )
! 953: return [D,M];
! 954: D *= CG; M *= CG;
! 955: break;
! 956: }
! 957: }
! 958: if ( U )
! 959: G = U;
! 960: else {
! 961: D += dp_hm(G); G = dp_rest(G);
! 962: }
! 963: }
! 964: return [D,M];
! 965: }
! 966:
! 967: def weyl_nf_quo_check(G,PS,R)
! 968: {
! 969: D = R[0]; M = R[1]; Coef = R[2];
! 970: Len = length(PS);
! 971: T = 0;
! 972: for ( I = 0; I < Len; I++ )
! 973: T += dp_weyl_mul(Coef[I],PS[I]);
! 974: return (M*G-T)==D;
! 975: }
! 976:
! 977: def weyl_nf_quo(B,G,M,PS)
! 978: {
! 979: Len = length(PS);
! 980: Coef = vector(Len);
! 981: for ( D = 0; G; ) {
! 982: for ( U = 0, L = B; L != []; L = cdr(L) ) {
! 983: if ( dp_redble(G,R=PS[car(L)]) > 0 ) {
! 984: GCD = igcd(dp_hc(G),dp_hc(R));
! 985: CG = idiv(dp_hc(R),GCD); CR = idiv(dp_hc(G),GCD);
! 986: for ( I = 0; I < Len; I++ ) Coef[I] *= CG;
! 987: Q = CR*dp_subd(G,R);
! 988: Coef[car(L)] += Q;
! 989: U = CG*G-dp_weyl_mul(Q,R);
! 990: D *= CG; M *= CG;
! 991: if ( !U )
! 992: return [D,M,Coef];
! 993: break;
! 994: }
! 995: }
! 996: if ( U )
! 997: G = U;
! 998: else {
! 999: D += dp_hm(G); G = dp_rest(G);
! 1000: }
! 1001: }
! 1002: return [D,M,Coef];
! 1003: }
! 1004:
! 1005: def weyl_nf_mod(B,G,PS,Mod)
! 1006: {
! 1007: for ( D = 0; G; ) {
! 1008: for ( U = 0, L = B; L != []; L = cdr(L) ) {
! 1009: if ( dp_redble(G,R=PS[car(L)]) > 0 ) {
! 1010: CR = dp_hc(G)/dp_hc(R);
! 1011: U = G-dp_weyl_mul_mod(CR*dp_mod(dp_subd(G,R),Mod,[]),R,Mod);
! 1012: if ( !U )
! 1013: return D;
! 1014: break;
! 1015: }
! 1016: }
! 1017: if ( U )
! 1018: G = U;
! 1019: else {
! 1020: D += dp_hm(G); G = dp_rest(G);
! 1021: }
! 1022: }
! 1023: return D;
! 1024: }
! 1025:
! 1026: def b_subst(F,V)
! 1027: {
! 1028: D = deg(F,V);
! 1029: C = newvect(D+1);
! 1030: for ( I = D; I >= 0; I-- )
! 1031: C[I] = coef(F,I,V);
! 1032: for ( I = 0, R = 0; I <= D; I++ )
! 1033: if ( C[I] )
! 1034: R += C[I]*v_factorial(V,I);
! 1035: return R;
! 1036: }
! 1037:
! 1038: def v_factorial(V,N)
! 1039: {
! 1040: for ( J = N-1, R = 1; J >= 0; J-- )
! 1041: R *= V-J;
! 1042: return R;
! 1043: }
! 1044:
! 1045: def w_tdeg(F,V,W)
! 1046: {
! 1047: dp_set_weight(newvect(length(W),W));
! 1048: T = dp_ptod(F,V);
! 1049: for ( R = 0; T; T = cdr(T) ) {
! 1050: D = dp_td(T);
! 1051: if ( D > R ) R = D;
! 1052: }
! 1053: return R;
! 1054: }
! 1055:
! 1056: def replace_vars_f(F)
! 1057: {
! 1058: return subst(F,s,TMP_S,t,TMP_T,y1,TMP_Y1,y2,TMP_Y2);
! 1059: }
! 1060:
! 1061: def replace_vars_v(V)
! 1062: {
! 1063: V = replace_var(V,s,TMP_S);
! 1064: V = replace_var(V,t,TMP_T);
! 1065: V = replace_var(V,y1,TMP_Y1);
! 1066: V = replace_var(V,y2,TMP_Y2);
! 1067: return V;
! 1068: }
! 1069:
! 1070: def replace_var(V,X,Y)
! 1071: {
! 1072: V = reverse(V);
! 1073: for ( R = []; V != []; V = cdr(V) ) {
! 1074: Z = car(V);
! 1075: if ( Z == X )
! 1076: R = cons(Y,R);
! 1077: else
! 1078: R = cons(Z,R);
! 1079: }
! 1080: return R;
! 1081: }
! 1082:
! 1083: def action_on_gfs(P,V,GFS)
! 1084: {
! 1085: DP = dp_ptod(P,V);
! 1086: N = length(V)/2;
! 1087: for ( I = N-1, V0 = []; I >= 0; I-- )
! 1088: V0 = cons(V[I],V0);
! 1089: R = [];
! 1090: for ( T = DP; T; T = dp_rest(T) )
! 1091: R = cons(action_on_gfs_1(dp_hm(T),N,V0,GFS),R);
! 1092: D = coef(car(R)[2],0);
! 1093: for ( T = cdr(R); T != []; T = cdr(T) ) {
! 1094: Di = coef(car(T)[2],0);
! 1095: if ( Di < D )
! 1096: D = Di;
! 1097: }
! 1098: F = GFS[1];
! 1099: for ( T = R, G = 0; T != []; T = cdr(T) )
! 1100: G += car(T)[0]*F^(car(T)[2]-(s+D));
! 1101: while ( 1 ) {
! 1102: G1 = tdiv(G,F);
! 1103: if ( G1 ) {
! 1104: G = G1;
! 1105: D++;
! 1106: } else
! 1107: return [G,F,s+D];
! 1108: }
! 1109: }
! 1110:
! 1111: def action_on_gfs_1(M,N,V,GFS)
! 1112: {
! 1113: G = GFS[0];
! 1114: F = GFS[1];
! 1115: S = GFS[2];
! 1116: C = dp_hc(M);
! 1117: E = dp_etov(M);
! 1118: for ( I = 0; I < N; I++ ) {
! 1119: VI = V[I];
! 1120: C *= VI^E[I];
! 1121: DFVI = diff(F,VI);
! 1122: for ( J = 0, EI = E[I+N]; J < EI; J++, S-- )
! 1123: G = diff(G,VI)*F+S*G*DFVI;
! 1124: }
! 1125: return [C*G,F,S];
! 1126: }
! 1127:
! 1128: /* stratification */
! 1129:
! 1130: def weyl_subst(F,P,V)
! 1131: {
! 1132: VF = var(F);
! 1133: D = deg(F,VF);
! 1134: P = dp_ptod(P,V);
! 1135: One = dp_ptod(1,V);
! 1136: for ( R = 0, I = D; I >= 0; I-- )
! 1137: R = dp_weyl_mul(R,P)+coef(F,I,VF)*One;
! 1138: return dp_dtop(R,V);
! 1139: }
! 1140:
! 1141: def bfactor(F)
! 1142: {
! 1143: L=length(F);
! 1144: for(I=0,B=1;I<L;I++)B*=F[I][0]^F[I][1];
! 1145: return fctr(B);
! 1146: }
! 1147:
! 1148: def gen_a(K)
! 1149: {
! 1150: D = x^(K+1);
! 1151: W = [];
! 1152: for ( I = 1; I <= K; I++ ) {
! 1153: D += (V=strtov("u"+rtostr(K-I+1)))*x^(K-I);
! 1154: W = cons(I+1,cons(V,W));
! 1155: }
! 1156: F = res(x,D,diff(D,x));
! 1157: return [D,F,reverse(W)];
! 1158: }
! 1159:
! 1160: def gen_d(K)
! 1161: {
! 1162: D = x^2*y+y^(K-1)+u1+u2*x+u3*x^2;
! 1163: W = reverse([u1,2*K-2,u2,K,u3,2]);
! 1164: U = [u3,u2,u1];
! 1165: for ( I = 4; I <= K; I++ ) {
! 1166: D += (V=strtov("u"+rtostr(I)))*y^(I-3);
! 1167: W = cons((2*K-2)-2*(I-3),cons(V,W));
! 1168: U = cons(V,U);
! 1169: }
! 1170: B = [D,diff(D,x),diff(D,y)];
! 1171: G = nd_gr_trace(B,append([x,y],U),1,1,0);
! 1172: G = nd_gr_trace(G,append([x,y],U),1,-1,[[0,2],[0,K]]);
! 1173: E = elimination(G,U);
! 1174: F = E[0];
! 1175: return [D,F,reverse(W)];
! 1176: }
! 1177:
! 1178: def gen_dm(K)
! 1179: {
! 1180: D = x^2*y-y^(K-1)+u1+u2*x+u3*x^2;
! 1181: W = reverse([u1,2*K-2,u2,K,u3,2]);
! 1182: U = [u3,u2,u1];
! 1183: for ( I = 4; I <= K; I++ ) {
! 1184: D += (V=strtov("u"+rtostr(I)))*y^(I-3);
! 1185: W = cons((2*K-2)-2*(I-3),cons(V,W));
! 1186: U = cons(V,U);
! 1187: }
! 1188: B = [D,diff(D,x),diff(D,y)];
! 1189: G = nd_gr_trace(B,append([x,y],U),1,1,0);
! 1190: G = nd_gr_trace(G,append([x,y],U),1,-1,[[0,2],[0,K]]);
! 1191: E = elimination(G,U);
! 1192: F = E[0];
! 1193: return [D,F,reverse(W)];
! 1194: }
! 1195:
! 1196: def elimination(G,V)
! 1197: {
! 1198: ANS=[];
! 1199: NG=length(G);
! 1200:
! 1201: for (I=NG-1;I>=0;I--)
! 1202: {
! 1203: VSet=vars(G[I]);
! 1204: DIFF=setminus(VSet,V);
! 1205:
! 1206: if ( DIFF ==[] )
! 1207: {
! 1208: ANS=cons(G[I],ANS);
! 1209: }
! 1210: }
! 1211: return ANS;
! 1212: }
! 1213:
! 1214: def weyl_ideal_quotient(B,F,VDV)
! 1215: {
! 1216: T = ttttt; DT = dttttt;
! 1217: J = cons((1-T)*F,vtol(ltov(B)*T));
! 1218: N = length(VDV); N1 = N/2;
! 1219: for ( I = N1-1, V1 = []; I >= 0; I-- )
! 1220: V1 = cons(VDV[I],V1);
! 1221: for ( I = 0, VDV1 = VDV; I < N1; I++ ) VDV1 = cdr(VDV1);
! 1222: VDV1 = append(cons(T,V1),cons(DT,VDV1));
! 1223: O1 = [[0,1],[0,N+1]];
! 1224: GJ = nd_weyl_gr(J,VDV1,0,O1);
! 1225: R = elimination(GJ,VDV);
! 1226: return R;
! 1227: }
! 1228:
! 1229: def bf_strat(F)
! 1230: {
! 1231: dp_ord(0);
! 1232: T0 = time();
! 1233: if ( type(Heu=getopt(heuristic)) == -1 ) Heu = 0;
! 1234: if ( type(Vord=getopt(vord)) == -1 || type(Vord) != 4 ) Vord = 0;
! 1235: if ( type(Wt=getopt(weight)) == -1 ) Wt = 0;
! 1236: L = in_ww(F|weight=Wt,heuristic=Heu,vord=Vord);
! 1237: T1 = time();
! 1238: print(["in_ww",(T1[0]+T1[1])-(T0[0]+T0[1])]);
! 1239:
! 1240: /* shortcuts */
! 1241: V = vars(F);
! 1242: dp_set_weight(0);
! 1243: dp_ord(0);
! 1244: Sing = sing(F,V);
! 1245: if ( Sing[0] == 1 || Sing[0] == -1 ) {
! 1246: return [[[F],[1],[[s+1,1]]],[[0],[F],[]]];
! 1247: } else if ( zero_dim(Sing,V,0) ) {
! 1248: N = length(V);
! 1249: P0 = [];
! 1250: for ( I = 0; I < N; I++ ) {
! 1251: M = minipoly(Sing,V,0,V[I],TMP_S);
! 1252: MF = cdr(fctr(M));
! 1253: if ( length(MF) == 1 && deg(MF[0][0],TMP_S)==1 ) {
! 1254: P0 = cons(subst(MF[0][0],TMP_S,V[I]),P0);
! 1255: } else break;
! 1256: }
! 1257: if ( I == N ) {
! 1258: Indata = L[0]; AllData = L[1]; VData = L[2];
! 1259: GIN = Indata[0]; VDV = Indata[1]; WVDV = AllData[4];
! 1260: W = Indata[4];
! 1261: dp_set_weight(W);
! 1262: B = weyl_minipoly(GIN,VDV,0,WVDV);
! 1263: B = subst(B,s,-s-1);
! 1264: dp_set_weight(0);
! 1265: return [
! 1266: [P0,[1],cdr(fctr(B))],
! 1267: [[F],P0,[[s+1,1]]],
! 1268: [[0],[F],[]]
! 1269: ];
! 1270: }
! 1271: }
! 1272:
! 1273: L2 = bf_strat_stage2(L);
! 1274: S = bf_strat_stage3(L2);
! 1275: R = [];
! 1276: for ( T = S; T != []; T = cdr(T) ) {
! 1277: Str = car(T);
! 1278: B = Str[2];
! 1279: B1 = [];
! 1280: for ( U = B; U != []; U = cdr(U) )
! 1281: B1 = cons([subst(car(U)[0],s,-s-1),car(U)[1]],B1);
! 1282: B1 = reverse(B1);
! 1283: R = cons([Str[0],Str[1],B1],R);
! 1284: }
! 1285: return reverse(R);
! 1286: }
! 1287:
! 1288: /* returns [QQ,V,V0,B,BF,W] */
! 1289: /* QQ : ideal in C[x,s] (s=tdt), V=[x1,..,xn,t], V0 = [x1,..,xn] */
! 1290: /* B : global b-function, BF : factor list of B, W : weight */
! 1291:
! 1292: def bf_strat_stage2(L)
! 1293: {
! 1294: T0 = time();
! 1295: InData = L[0]; VData = L[2];
! 1296: G1 = InData[0]; VDV = InData[1]; W = InData[4]; W0 = VData[4];
! 1297: N = length(VDV); N1 = N/2;
! 1298: V = InData[2]; DV = InData[3];
! 1299: T = VData[2]; DT = VData[3];
! 1300: V0 = VData[0]; DVR = VData[1];
! 1301: dp_set_weight(W);
! 1302: for ( I = 0; DVR != []; I++, DVR = cdr(DVR) ) {
! 1303: DVRV = cons(DT,append(cdr(DVR),V));
! 1304: M = elim_mat(VDV,DVRV);
! 1305: for ( K = 0; K < N; K++ )
! 1306: M[1][K] = W[K];
! 1307: dp_ord(0); H1 = map(dp_ht,map(dp_ptod,G1,VDV));
! 1308: dp_ord(M); H2 = map(dp_ht,map(dp_ptod,G1,VDV));
! 1309: if ( H1 == H2 )
! 1310: G2 = G1;
! 1311: else
! 1312: G2 = nd_weyl_gr_trace(G1,VDV,0,-1,M);
! 1313: G1 = elimination(G2,DVRV);
! 1314: }
! 1315: T1 = time();
! 1316: B = weyl_minipoly(G1,VDV,0,T*DT);
! 1317: T2 = time();
! 1318: BF = cdr(fctr(B));
! 1319:
! 1320: dp_set_weight(0);
! 1321: G1 = map(psi0,G1,T,DT);
! 1322: QQ = map(subst,map(b_subst,map(subst,G1,DT,1),T),T,var(B));
! 1323: if ( type(getopt(ideal)) != -1 ) return [QQ,V];
! 1324: print(["elim",(T1[0]+T1[1])-(T0[0]+T0[1])]);
! 1325: print(["globalb",(T2[0]+T2[1])-(T1[0]+T1[1])]);
! 1326: return [QQ,V,V0,B,BF,W0,DV];
! 1327: }
! 1328:
! 1329: def bf_strat_stage3(L)
! 1330: {
! 1331: T0 = time();
! 1332: QQ = L[0]; V0 = L[2]; B = L[3]; BF = L[4]; W0 = L[5];
! 1333: NF = length(BF);
! 1334: Data = vector(NF);
! 1335: for ( I = J = 0; I < NF; I++ ) {
! 1336: DI = tower_in_p(QQ,B,BF[I],V0,W0);
! 1337: NDI = length(DI);
! 1338: for ( K = 0; K < J; K++ ) {
! 1339: if ( length(DK=Data[K]) == NDI ) {
! 1340: for ( L = 0; L < NDI; L++ ) {
! 1341: CL = DI[L][1]; CK = DK[L][1];
! 1342: if ( !zero_inclusion(CL,CK,V0)
! 1343: || !zero_inclusion(CK,CL,V0) ) break;
! 1344: }
! 1345: if ( L == NDI ) break;
! 1346: }
! 1347: }
! 1348: if ( K < J ) {
! 1349: for ( L = 0, T = []; L < NDI; L++ )
! 1350: T = cons([[DK[L][0][0]*DI[L][0][0],DK[L][0][1]],
! 1351: DK[L][1],DK[L][2]],T);
! 1352: Data[K] = reverse(T);
! 1353: } else
! 1354: Data[J++] = DI;
! 1355: }
! 1356: Data1 = vector(J);
! 1357: for ( I = 0; I < J; I++ )
! 1358: Data1[I] = Data[I];
! 1359: T1 = time();
! 1360: Str = stratify_bf(Data1,V0);
! 1361: T2 = time();
! 1362: print(["tower",(T1[0]+T1[1])-(T0[0]+T0[1])]);
! 1363: print(["strat",(T2[0]+T2[1])-(T1[0]+T1[1])]);
! 1364: return Str;
! 1365: }
! 1366:
! 1367: /*
! 1368: InData = [GIN,VDV,V,DV,WtV]
! 1369: AllData = [G0,GIN0,VDV0,W,WVDV,WtV0]
! 1370: */
! 1371:
! 1372: def bf_local(F,P)
! 1373: {
! 1374: if ( type(Heu=getopt(heuristic)) == -1 ) Heu = 0;
! 1375: if ( type(Vord=getopt(vord)) == -1 || type(Vord) != 4 ) Vord = 0;
! 1376: if ( type(Wt=getopt(weight)) == -1 ) Wt = 0;
! 1377: L = in_ww(F|weight=Wt,heuristic=Heu,vord=Vord);
! 1378: InData = L[0]; AllData = L[1]; VData = L[2];
! 1379: G = InData[0]; VDV = InData[1];
! 1380: V = InData[2]; DV = InData[3];
! 1381:
! 1382: V0 = VData[0]; DV0 = VData[1]; T = VData[2]; DT = VData[3]; W0 = VData[4];
! 1383:
! 1384: L2 = bf_strat_stage2(L);
! 1385:
! 1386: /* L2 = [QQ,V,V0,B,BF,W] */
! 1387: /* QQ : ideal in C[x,s] (s=tdt), V=[t,x1,..,xn], V0 = [x1,..,xn] */
! 1388: /* B : global b-function, BF : factor list of B, W : weight */
! 1389:
! 1390: QQ = L2[0]; B = L2[3]; BF = L2[4]; W = L2[5];
! 1391:
! 1392: NF = length(BF);
! 1393: BP = [];
! 1394: dp_set_weight(0);
! 1395: for ( I = J = 0; I < NF; I++ ) {
! 1396: List = compute_exponent(QQ,B,BF[I],P,V0,W0);
! 1397: DI = List[0]; QQI = List[1];
! 1398: if ( DI )
! 1399: BP = cons([BF[I][0],DI],BP);
! 1400: if ( I == 0 )
! 1401: Id = QQI;
! 1402: else
! 1403: Id = ideal_intersection(Id,QQI,V0,0);
! 1404: }
! 1405: for ( List = Id; List != []; List = cdr(List) )
! 1406: if ( subst_vars(car(List),P) )
! 1407: break;
! 1408: if ( List == [] ) error("bf_local : inconsitent intersection");
! 1409: Ax = car(List);
! 1410: for ( BPT = 1, List = BP; List != []; List = cdr(List) )
! 1411: BPT *= car(List)[0]^car(List)[1];
! 1412: BPT = weyl_subst(BPT,T*DT,VDV);
! 1413:
! 1414: /* computation using G0,GIN0,VDV0 */
! 1415: G0 = AllData[0]; GIN0 = AllData[1]; VDV0 = AllData[2]; WtV0 = AllData[5];
! 1416: dp_set_weight(WtV0); dp_ord(0);
! 1417: PS = map(dp_ptod,GIN0,VDV0); Len = length(PS);
! 1418: for ( I = Len-1, Ind = []; I >= 0; I-- ) Ind = cons(I,Ind);
! 1419: /* QR = [D,M,Coef] */
! 1420: AxBPT = dp_ptod(Ax*BPT,VDV0);
! 1421: QR = weyl_nf_quo(Ind,AxBPT,1,PS);
! 1422: if ( !weyl_nf_quo_check(AxBPT,PS,QR) ) error("bf_local : invalid quotient");
! 1423: if ( QR[0] ) error("bf_local : invalid quotient");
! 1424: Den = QR[1]; Coef = QR[2];
! 1425: for ( I = 0, R = Den*AxBPT; I < Len; I++ )
! 1426: R -= dp_weyl_mul(Coef[I],dp_ptod(G0[I],VDV0));
! 1427: R = dp_dtop(R,VDV0);
! 1428: CR = conv_tdt(R,F,V0,DV0,T,DT);
! 1429:
! 1430: dp_set_weight(0);
! 1431: return [BP,Ax,CR];
! 1432: }
! 1433:
! 1434: /* t^(l+k)*dt^l (k>l) -> (s-k)(s-k-1)...(s-(k+l-1))t^k */
! 1435: def conv_tdt(P,F,V0,DV0,T,DT)
! 1436: {
! 1437: DP = dp_ptod(P,[T,DT]);
! 1438: VDV = append(cons(T,V0),cons(DT,DV0));
! 1439: R = 0;
! 1440: DF = dp_ptod(F,VDV);
! 1441: for ( ; DP; DP = dp_rest(DP) ) {
! 1442: C = dp_hc(DP);
! 1443: E = dp_etov(dp_ht(DP));
! 1444: L = E[1]; K = E[0]-E[1];
! 1445: S = 1;
! 1446: for ( I = 0; I < L; I++ )
! 1447: S *= ((-T-1)-K-I);
! 1448: U = dp_ptod(C*S,VDV);
! 1449: for ( I = 1; I < K; I++ )
! 1450: U = dp_weyl_mul(U,DF);
! 1451: R += dp_dtop(U,VDV);
! 1452: }
! 1453: return subst(R,T,s);
! 1454: }
! 1455:
! 1456: def merge_tower(Str,Tower,V)
! 1457: {
! 1458: Prev = car(Tower); T = cdr(Tower);
! 1459: Str1 = [];
! 1460: for ( ; T != []; T = cdr(T) ) {
! 1461: Cur = car(T);
! 1462: Str1 = cons([Cur[1],Prev[1],[Prev[0]]],Str1);
! 1463: Prev = Cur;
! 1464: }
! 1465: Str1 = cons([[0],Prev[1],[]],Str1);
! 1466: Str1 = reverse(Str1);
! 1467: if ( Str == [] ) return Str1;
! 1468: U = [];
! 1469: for ( T = Str; T != []; T = cdr(T) ) {
! 1470: for ( S = Str1; S != []; S = cdr(S) ) {
! 1471: Int = int_cons(T[0],S[0],V);
! 1472: if ( Int[0] != [1] )
! 1473: U = cons(append(Int,[append(T[0][2],S[0][2])]),U);
! 1474: }
! 1475: }
! 1476: U = reverse(U);
! 1477: return U;
! 1478: }
! 1479:
! 1480: def stratify_bf(Data,V)
! 1481: {
! 1482: N = length(Data);
! 1483: R = [];
! 1484: for ( I = 0; I < N; I++ )
! 1485: R = merge_tower(R,Data[I],V);
! 1486: return R;
! 1487: }
! 1488:
! 1489: def tdt_homogenize(F,L)
! 1490: {
! 1491: TY1 = L[0]; T = TY1[0]; Y1 = TY1[1];
! 1492: TY2 = L[1]; DT = TY2[0]; Y2 = TY2[1];
! 1493: DF = dp_ptod(F,[T,DT]);
! 1494: for ( R = 0; DF; DF = dp_rest(DF) ) {
! 1495: M = dp_hm(DF);
! 1496: E = dp_etov(M);
! 1497: W = E[1]-E[0];
! 1498: if ( W > 0 ) R += Y1^W*dp_dtop(M,[T,DT]);
! 1499: else R += Y2^W*dp_dtop(M,[T,DT]);
! 1500: }
! 1501: return R;
! 1502: }
! 1503:
! 1504: def sing(F,V)
! 1505: {
! 1506: R = [F];
! 1507: for ( T = V; T != []; T = cdr(T) )
! 1508: R = cons(diff(F,car(T)),R);
! 1509: G = nd_gr_trace(R,V,1,1,0);
! 1510: return G;
! 1511: }
! 1512:
! 1513: def tower_in_p(B,F,FD,V,W)
! 1514: {
! 1515: TT = ttttt;
! 1516: N = length(V); S = var(F); SV = cons(S,V); V1 = cons(TT,SV);
! 1517: Wt = append(append([1,1],W),[1]);
! 1518: dp_set_weight(Wt);
! 1519:
! 1520: F1 = FD[0]; D = FD[1];
! 1521: O1 = [[0,1],[0,N+1]]; O2 = [[0,1],[0,N]];
! 1522:
! 1523: TF = sdiv(F,F1^D);
! 1524:
! 1525: T = nd_gr_trace(cons((1-TT)*TF,vtol(TT*ltov(B))),V1,1,1,O1);
! 1526: T = elimination(T,SV);
! 1527: Q = map(sdiv,T,TF);
! 1528: dp_set_weight(cdr(Wt));
! 1529: QQ = elimination(nd_gr(Q,SV,0,O2),V);
! 1530: E = [[[F1,0],QQ,PD]];
! 1531:
! 1532: for ( I = D-1; I >= 0; I-- ) {
! 1533: dp_set_weight(Wt);
! 1534: T = nd_gr_trace(cons((1-TT)*F1,vtol(TT*ltov(Q))),V1,1,1,O1);
! 1535: T = elimination(T,SV);
! 1536: Q = map(sdiv,T,F1);
! 1537: dp_set_weight(cdr(Wt));
! 1538: QQ = elimination(nd_gr(Q,SV,0,O2),V);
! 1539: E = cons([[F1,D-I],QQ,PD],E);
! 1540: }
! 1541: dp_set_weight(0);
! 1542: return E;
! 1543: }
! 1544:
! 1545: def subst_vars(F,P)
! 1546: {
! 1547: for ( ; P != []; P = cdr(cdr(P)) )
! 1548: F = subst(F,P[0],P[1]);
! 1549: return F;
! 1550: }
! 1551:
! 1552: def compute_exponent(B,F,FD,P,V,W)
! 1553: {
! 1554: TT = ttttt;
! 1555: N = length(V); S = var(F); SV = cons(S,V); V1 = cons(TT,SV);
! 1556: F1 = FD[0]; D = FD[1];
! 1557:
! 1558: Wt = append(append([1,1],W),[1]);
! 1559: dp_set_weight(Wt);
! 1560: O1 = [[0,1],[0,N+1]]; O2 = [[0,1],[0,N]];
! 1561:
! 1562: TF = sdiv(F,F1^D);
! 1563:
! 1564: T = nd_gr_trace(cons((1-TT)*TF,vtol(TT*ltov(B))),V1,0,1,O1);
! 1565: T = elimination(T,SV);
! 1566: Q = map(sdiv,T,TF);
! 1567: QQ = elimination(nd_gr(Q,SV,0,O2),V);
! 1568:
! 1569: for ( I = 0; I < D; I++ ) {
! 1570: for ( T = QQ; T != []; T = cdr(T) )
! 1571: if ( subst_vars(car(T),P) ) {
! 1572: dp_set_weight(0);
! 1573: return [I,QQ];
! 1574: }
! 1575: T = nd_gr_trace(cons((1-TT)*F1,vtol(TT*ltov(Q))),V1,0,1,O1);
! 1576: T = elimination(T,SV);
! 1577: Q = map(sdiv,T,F1);
! 1578: QQ = elimination(nd_gr(Q,SV,0,O2),V);
! 1579: }
! 1580: dp_set_weight(0);
! 1581: return [D,QQ];
! 1582: }
! 1583:
! 1584: /* V(B) subset V(A) ? */
! 1585:
! 1586: def zero_inclusion(A,B,V)
! 1587: {
! 1588: NV = ttttt;
! 1589: for ( T = A; T != []; T = cdr(T) ) {
! 1590: F = car(T);
! 1591: G = nd_gr_trace(cons(NV*F-1,B),cons(NV,V),1,-1,0);
! 1592: if ( type(car(G)) != 1 ) return 0;
! 1593: }
! 1594: return 1;
! 1595: }
! 1596:
! 1597: def weyl_divide_by_right(G,H,V,O)
! 1598: {
! 1599: dp_ord(O); G = dp_ptod(G,V); H = dp_ptod(H,V);
! 1600: CH = dp_hc(H);
! 1601: for ( Q = 0; G; ) {
! 1602: if ( !dp_redble(G,H) ) return 0;
! 1603: CG = dp_hc(G);
! 1604: Q1 = CG/CH*dp_subd(G,H);
! 1605: G -= dp_weyl_mul(Q1,H);
! 1606: Q += Q1;
! 1607: }
! 1608: return dp_dtop(Q,V);
! 1609: }
! 1610:
! 1611: def elim_mat(V,W)
! 1612: {
! 1613: N = length(V);
! 1614: M = matrix(N+1,N);
! 1615: for ( J = 0; J < N; J++ ) if ( !member(V[J],W) ) M[0][J] = 1;
! 1616: for ( J = 0; J < N; J++ ) M[1][J] = 1;
! 1617: for ( I = 2; I <= N; I++ ) M[I][N-I+1] = -1;
! 1618: return M;
! 1619: }
! 1620:
! 1621: /* (P-Q)cap(R-S)=(P cap Q^c)cap(R cap S^c)=(P cap R)cap(Q cup S)^c
! 1622: =(P cap R)-(Q cup S)
! 1623: */
! 1624:
! 1625: def int_cons(A,B,V)
! 1626: {
! 1627: AZ = A[0]; AN = A[1];
! 1628: BZ = B[0]; BN = B[1];
! 1629: CZ = nd_gr_trace(append(AZ,BZ),V,1,1,0);
! 1630: CN = ideal_intersection(AN,BN,V,0);
! 1631: if ( zero_inclusion(CN,CZ,V) )
! 1632: return [[1],[]];
! 1633: else
! 1634: return [CZ,CN];
! 1635: }
! 1636:
! 1637: def ideal_intersection(A,B,V,Ord)
! 1638: {
! 1639: T = ttttt;
! 1640: G = nd_gr_trace(append(vtol(ltov(A)*T),vtol(ltov(B)*(1-T))),
! 1641: cons(T,V),1,1,[[0,1],[Ord,length(V)]]);
! 1642: return elimination(G,V);
! 1643: }
! 1644: endmodule $
! 1645: end$
! 1646:
! 1647:
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>
![[BACK]](/icons/cvsweb/back.gif){kind=icon}
Return to ndbf.rr CVS log ![[TXT]](/icons/cvsweb/text.gif){kind=icon}
![[DIR]](/icons/cvsweb/dir.gif){kind=icon}
Up to [local] / OpenXM / src / asir-contrib / testing / noro