![[BACK]](/icons/cvsweb/back.gif){kind=icon}
Return to pd.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|
Return to pd.rr CVS log | Up to [local] / OpenXM / src / asir-contrib / testing / noro |
1.1 noro 1: import("gr")$
2: module noro_pd$
3:
4: static GBCheck,F4,Procs,SatHomo$
5:
6: localf init_procs, kill_procs, syca_dec, syc_dec, find_separating_ideal0$
7: localf find_separating_ideal1, find_separating_ideal2$
8: localf sy_dec, pseudo_dec, iso_comp, prima_dec$
9: localf prime_dec, prime_dec_main, lex_predec1, zprimedec, zprimadec$
10: localf complete_qdecomp, partial_qdecomp, partial_qdecomp0, complete_decomp$
11: localf partial_decomp, partial_decomp0, zprimacomp, zprimecomp$
1.5 ! noro 12: localf fast_gb, incremental_gb, elim_gb, ldim, make_mod_subst$
1.1 noro 13: localf rsgn, find_npos, gen_minipoly, indepset$
14: localf maxindep, contraction, ideal_list_intersection, ideal_intersection$
15: localf radical_membership, quick_radical_membership, modular_radical_membership$
16: localf radical_membership_rep, ideal_product, saturation$
17: localf sat, satind, sat_ind, colon$
1.2 noro 18: localf ideal_colon, ideal_sat, ideal_inclusion, qd_simp_comp, qd_remove_redundant_comp$
1.1 noro 19: localf remove_redundant_comp, remove_redundant_comp_first, ppart, sq$
20: localf lcfactor, compute_deg0, compute_deg, member$
21: localf elimination, setintersection, setminus, sep_list$
22: localf first_element, comp_tdeg, tdeg, comp_by_ord, comp_by_second$
23: localf gbcheck,f4,sathomo$
24:
25: SatHomo=0$
26: GBCheck=1$
27:
28: #define MAX(a,b) ((a)>(b)?(a):(b))
29:
30: def gbcheck(A)
31: {
32: if ( A ) GBCheck = 1;
1.3 noro 33: else GBCheck = -1;
1.1 noro 34: }
35:
36: def f4(A)
37: {
38: if ( A ) F4 = 1;
39: else F4 = 0;
40: }
41:
42: def sathomo(A)
43: {
44: if ( A ) SatHomo = 1;
45: else SatHomo = 0;
46: }
47:
48: def init_procs()
49: {
50: if ( type(NoX=getopt(nox)) == -1 ) NoX = 0;
51: if ( !Procs ) {
52: if ( NoX ) {
53: P0 = ox_launch_nox();
54: P1 = ox_launch_nox();
55: } else {
56: P0 = ox_launch();
57: P1 = ox_launch();
58: }
59: Procs = [P0,P1];
60: }
61: }
62:
63: def kill_procs()
64: {
65: if ( Procs ) {
66: ox_shutdown(Procs[0]);
67: ox_shutdown(Procs[1]);
68: Procs = 0;
69: }
70: }
71:
72: /* SYC primary decomositions */
73:
74: def syca_dec(B,V)
75: {
76: T00 = time();
77: if ( type(Nolexdec=getopt(nolexdec)) == -1 ) Nolexdec = 0;
78: if ( type(SepIdeal=getopt(sepideal)) == -1 ) SepIdeal = 1;
79: if ( type(NoSimp=getopt(nosimp)) == -1 ) NoSimp = 0;
80: if ( type(Time=getopt(time)) == -1 ) Time = 0;
81: Ord = 0;
82: Gt = G0 = G = fast_gb(B,V,0,Ord);
83: Q0 = Q = []; IntQ0 = IntQ = [1]; First = 1;
84: C = 0;
85:
86: Tass = Tiso = Tcolon = Tsep = Tirred = 0;
87: Rass = Riso = Rcolon = Rsep = Rirred = 0;
88: while ( 1 ) {
89: if ( type(Gt[0])==1 ) break;
90: T0 = time();
91: Pt = prime_dec(Gt,V|indep=1,nolexdec=Nolexdec);
92: T1 = time(); Tass += T1[0]-T0[0]+T1[1]-T0[1]; Rass += T1[3]-T0[3];
93: T0 = time();
94: Qt = iso_comp(Gt,Pt,V,Ord);
95: T1 = time(); Tiso += T1[0]-T0[0]+T1[1]-T0[1]; Riso += T1[3]-T0[3];
96: IntQt = ideal_list_intersection(map(first_element,Qt),V,Ord);
97: IntPt = ideal_list_intersection(map(first_element,Pt),V,Ord);
98: if ( First ) {
99: IntQ0 = IntQ = IntQt; IntP = IntPt; Qi = Qt; First = 0;
100: } else {
101: IntQ1 = ideal_intersection(IntQ,IntQt,V,Ord);
102: if ( gb_comp(IntQ,IntQ1) ) {
103: G = Gt; IntP = IntPt; Q = []; IntQ = [1]; C = 0;
104: continue;
105: } else {
106: IntQ = IntQ1;
107: IntQ1 = ideal_intersection(IntQ0,IntQt,V,Ord);
108: if ( !gb_comp(IntQ0,IntQ1) ) {
109: IntQ0 = IntQ1;
110: Q = append(Qt,Q); Q0 = append(Qt,Q0);
111: }
112: }
113: }
114: if ( gb_comp(IntQt,Gt) || gb_comp(IntQ,G) || gb_comp(IntQ0,G0) ) break;
115: T0 = time();
116: C1 = ideal_colon(G,IntQ,V);
117: T1 = time(); Tcolon += T1[0]-T0[0]+T1[1]-T0[1]; Rcolon += T1[3]-T0[3];
118: if ( C && gb_comp(C,C1) ) {
119: G = Gt; IntP = IntPt; Q = []; IntQ = [1]; C = 0;
120: continue;
121: } else C = C1;
122: T0 = time();
123: if ( SepIdeal == 0 )
124: Ok = find_separating_ideal0(C,G,IntQ,IntP,V,Ord);
125: else if ( SepIdeal == 1 )
126: Ok = find_separating_ideal1(C,G,IntQ,IntP,V,Ord);
127: else if ( SepIdeal == 2 )
128: Ok = find_separating_ideal2(C,G,IntQ,IntP,V,Ord);
129: G1 = append(Ok,G);
130: Gt1 = fast_gb(G1,V,0,Ord);
131: T1 = time(); Tsep += T1[0]-T0[0]+T1[1]-T0[1]; Rsep += T1[3]-T0[3];
132: #if 0
133: if ( ideal_inclusion(Gt1,Gt,V,Ord) ) {
134: G = Gt; IntP = IntPt; Q = []; IntQ = [1]; C = 0;
135: } else
136: #endif
137: Gt = Gt1;
138: }
139: T0 = time();
140: if ( !NoSimp ) Q1 = qd_remove_redundant_comp(G0,Qi,Q0,V,Ord);
141: else Q1 = Q0;
142: if ( Time ) {
143: T1 = time(); Tirred += T1[0]-T0[0]+T1[1]-T0[1]; Rirred += T1[3]-T0[3];
144: Tall = T1[0]-T00[0]+T1[1]-T00[1]; Rall += T1[3]-T00[3];
145: print(["total",Tall,"ass",Tass,"iso",Tiso, "colon",Tcolon,"sep",Tsep,"irred",Tirred]);
146: print(["Rtotal",Rall,"Rass",Rass,"Riso",Riso, "Rcolon",Rcolon,"Rsep",Rsep,"Rirred",Rirred]);
147: print(["iso",length(Qi),"emb",length(Q0),"->",length(Q1)]);
148: }
149: return [Qi,Q1];
150: }
151:
152: def syc_dec(B,V)
153: {
154: T00 = time();
155: if ( type(Nolexdec=getopt(nolexdec)) == -1 ) Nolexdec = 0;
156: if ( type(SepIdeal=getopt(sepideal)) == -1 ) SepIdeal = 1;
157: if ( type(NoSimp=getopt(nosimp)) == -1 ) NoSimp = 0;
158: if ( type(Time=getopt(time)) == -1 ) Time = 0;
159: Ord = 0;
160: G = fast_gb(B,V,0,Ord);
161: Q = []; IntQ = [1]; Gt = G; First = 1;
162: Tass = Tiso = Tcolon = Tsep = Tirred = 0;
163: Rass = Riso = Rcolon = Rsep = Rirred = 0;
164: while ( 1 ) {
165: if ( type(Gt[0])==1 ) break;
166: T0 = time();
167: Pt = prime_dec(Gt,V|indep=1,nolexdec=Nolexdec);
168: T1 = time(); Tass += T1[0]-T0[0]+T1[1]-T0[1]; Rass += T1[3]-T0[3];
169: T0 = time();
170: Qt = iso_comp(Gt,Pt,V,Ord);
171: T1 = time(); Tiso += T1[0]-T0[0]+T1[1]-T0[1]; Riso += T1[3]-T0[3];
172: IntQt = ideal_list_intersection(map(first_element,Qt),V,Ord);
173: IntPt = ideal_list_intersection(map(first_element,Pt),V,Ord);
174: if ( First ) {
175: IntQ = IntQt; Qi = Qt; First = 0;
176: } else {
177: IntQ1 = ideal_intersection(IntQ,IntQt,V,Ord);
178: if ( !gb_comp(IntQ1,IntQ) )
179: Q = append(Qt,Q);
180: }
181: if ( gb_comp(IntQ,G) || gb_comp(IntQt,Gt) )
182: break;
183: T0 = time();
184: C = ideal_colon(Gt,IntQt,V);
185: T1 = time(); Tcolon += T1[0]-T0[0]+T1[1]-T0[1]; Rcolon += T1[3]-T0[3];
1.2 noro 186: T0 = time();
1.1 noro 187: if ( SepIdeal == 0 )
188: Ok = find_separating_ideal0(C,Gt,IntQt,IntPt,V,Ord);
189: else if ( SepIdeal == 1 )
190: Ok = find_separating_ideal1(C,Gt,IntQt,IntPt,V,Ord);
191: else if ( SepIdeal == 2 )
192: Ok = find_separating_ideal2(C,Gt,IntQt,IntPt,V,Ord);
193: G1 = append(Ok,Gt);
194: Gt = fast_gb(G1,V,0,Ord);
195: T1 = time(); Tsep += T1[0]-T0[0]+T1[1]-T0[1]; Rsep += T1[3]-T0[3];
196: }
197: T0 = time();
198: if ( !NoSimp ) Q1 = qd_remove_redundant_comp(G,Qi,Q,V,Ord);
199: else Q1 = Q;
200: T1 = time(); Tirred += T1[0]-T0[0]+T1[1]-T0[1]; Rirred += T1[3]-T0[3];
201: Tall = T1[0]-T00[0]+T1[1]-T00[1]; Rall += T1[3]-T00[3];
202: if ( Time ) {
203: print(["total",Tall,"ass",Tass,"iso",Tiso, "colon",Tcolon,"sep",Tsep,"irred",Tirred]);
204: print(["Rtotal",Rall,"Rass",Rass,"Riso",Riso, "Rcolon",Rcolon,"Rsep",Rsep,"Rirred",Rirred]);
205: print(["iso",length(Qi),"emb",length(Q),"->",length(Q1)]);
206: }
207: return [Qi,Q1];
208: }
209:
210: /* C=G:Q, Rad=rad(Q), return J s.t. Q cap (G+J) = G */
211:
212: def find_separating_ideal0(C,G,Q,Rad,V,Ord) {
213: for ( CI = C, I = 1; ; I++ ) {
214: for ( T = CI, S = []; T != []; T = cdr(T) )
215: if ( nd_nf(car(T),Q,V,Ord,0) ) S = cons(car(T),S);
216: if ( S == [] )
217: error("find_separating_ideal0 : cannot happen");
218: G1 = append(S,G);
219: Int = ideal_intersection(G1,Q,V,Ord);
220: /* check whether (Q cap (G+S)) = G */
221: if ( gb_comp(Int,G) ) return reverse(S);
222: CI = ideal_product(CI,C,V);
223: }
224: }
225:
226: def find_separating_ideal1(C,G,Q,Rad,V,Ord) {
227: for ( T = C, S = []; T != []; T = cdr(T) )
228: if ( nd_nf(car(T),Q,V,Ord,0) ) S = cons(car(T),S);
229: if ( S == [] )
230: error("find_separating_ideal1 : cannot happen");
231: G1 = append(S,G);
232: Int = ideal_intersection(G1,Q,V,Ord);
233: /* check whether (Q cap (G+S)) = G */
234: if ( gb_comp(Int,G) ) return reverse(S);
235:
1.5 ! noro 236: /* or qsort(C,comp_tdeg) */
1.1 noro 237: C = qsort(S,comp_tdeg);
1.5 ! noro 238:
! 239: Tmp = ttttt; TV = cons(Tmp,V); Ord1 = [[0,1],[Ord,length(V)]];
! 240: Int0 = incremental_gb(append(vtol(ltov(G)*Tmp),vtol(ltov(Q)*(1-Tmp))),
! 241: TV,Ord1|gbblock=[[0,length(G)]]);
1.1 noro 242: for ( T = C, S = []; T != []; T = cdr(T) ) {
243: if ( !nd_nf(car(T),Rad,V,Ord,0) ) continue;
244: Ui = U = car(T);
245: for ( I = 1; ; I++ ) {
246: G1 = cons(Ui,G);
247: Int = ideal_intersection(G1,Q,V,Ord);
248: if ( gb_comp(Int,G) ) break;
249: else
250: Ui = nd_nf(Ui*U,G,V,Ord,0);
251: }
1.5 ! noro 252: Int1 = incremental_gb(append(Int0,[Tmp*Ui]),TV,Ord1
! 253: |gbblock=[[0,length(Int0)]]);
! 254: Int = elimination(Int1,V);
! 255: if ( !gb_comp(Int,G) )
! 256: break;
! 257: else {
! 258: Int0 = Int1;
! 259: S = cons(Ui,S);
1.1 noro 260: }
261: }
262: return reverse(S);
263: }
264:
265: def find_separating_ideal2(C,G,Q,Rad,V,Ord) {
266: for ( T = C, S = []; T != []; T = cdr(T) )
267: if ( nd_nf(car(T),Q,V,Ord,0) ) S = cons(car(T),S);
268: if ( S == [] )
269: error("find_separating_ideal2 : cannot happen");
270: G1 = append(S,G);
271: Int = ideal_intersection(G1,Q,V,Ord);
272: /* check whether (Q cap (G+S)) = G */
273: if ( gb_comp(Int,G) ) return reverse(S);
274:
1.5 ! noro 275: /* or qsort(S,comp_tdeg) */
1.1 noro 276: C = qsort(C,comp_tdeg);
1.5 ! noro 277: Dp = dp_gr_print(); dp_gr_print(0);
1.1 noro 278: for ( T = C, S = []; T != []; T = cdr(T) ) {
1.5 ! noro 279: print(length(T));
1.1 noro 280: if ( !nd_nf(car(T),Rad,V,Ord,0) ) continue;
281: Ui = U = car(T);
282: for ( I = 1; ; I++ ) {
283: G1 = cons(Ui,G);
284: Int = ideal_intersection(G1,Q,V,Ord);
285: if ( gb_comp(Int,G) ) break;
286: else
287: Ui = nd_nf(Ui*U,G,V,Ord,0);
288: }
289: S = cons(Ui,S);
290: }
1.5 ! noro 291: S = qsort(S,comp_tdeg);
! 292: /* S = reverse(S); */
1.1 noro 293: Len = length(S);
1.5 ! noro 294:
! 295: Tmp = ttttt; TV = cons(Tmp,V); Ord1 = [[0,1],[Ord,length(V)]];
1.1 noro 296: if ( Len > 1 ) {
1.5 ! noro 297: Prev = 1;
! 298: Cur = 2;
! 299: G1 = append(G,[S[0]]);
! 300: Int0 = incremental_gb(append(vtol(ltov(G1)*Tmp),vtol(ltov(Q)*(1-Tmp))),
! 301: TV,Ord1|gbblock=[[0,length(G)]]);
! 302: while ( Prev < Cur ) {
! 303: for ( St = [], I = Prev; I < Cur; I++ ) St = cons(Tmp*S[I],St);
! 304: Int1 = incremental_gb(append(Int0,St),TV,Ord1
! 305: |gbblock=[[0,length(Int0)]]);
! 306: Int = elimination(Int1,V);
! 307: if ( gb_comp(Int,G) ) {
! 308: print(Cur);
! 309: Prev = Cur;
! 310: Cur = Cur+idiv(Len-Cur+1,2);
! 311: Int0 = Int1;
! 312: } else {
! 313: Cur = Prev + idiv(Cur-Prev,2);
1.1 noro 314: }
315: }
1.5 ! noro 316: for ( St = [], I = 0; I < Prev; I++ ) St = cons(S[I],St);
! 317: Ok = reverse(St);
! 318: } else
! 319: Ok = [S[0]];
! 320: print([length(S),length(Ok)]);
! 321: dp_gr_print(Dp);
1.1 noro 322: return Ok;
323: }
324:
325: /* SY primary decompsition */
326:
327: def sy_dec(B,V)
328: {
329: if ( type(Nolexdec=getopt(nolexdec)) == -1 ) Nolexdec = 0;
330: Ord = 0;
331: G = fast_gb(B,V,0,Ord);
332: Q = [];
333: IntQ = [1];
334: Gt = G;
335: First = 1;
336: while ( 1 ) {
337: if ( type(Gt[0]) == 1 ) break;
338: Pt = prime_dec(Gt,V|indep=1,nolexdec=Nolexdec);
339: L = pseudo_dec(Gt,Pt,V,Ord);
340: Qt = L[0]; Rt = L[1]; St = L[2];
341: IntQt = ideal_list_intersection(Qt,V,Ord);
342: if ( First ) {
343: IntQ = IntQt;
344: Qi = Qt;
345: First = 0;
346: } else {
347: IntQ = ideal_intersection(IntQ,IntQt,V,Ord);
348: Q = append(Qt,Q);
349: }
350: if ( gb_comp(IntQ,G) ) break;
351: for ( T = Rt; T != []; T = cdr(T) ) {
352: if ( type(car(T)[0]) == 1 ) continue;
353: U = sy_dec(car(T),V|nolexdec=Nolexdec);
354: IntQ = ideal_list_intersection(cons(IntQ,U),V,Ord);
355: Q = append(U,Q);
356: if ( gb_comp(IntQ,G) ) break;
357: }
358: Gt = fast_gb(append(Gt,St),V,0,Ord);
359: }
360: Q = remove_redundant_comp(G,Qi,Q,V,Ord);
361: return append(Qi,Q);
362: }
363:
364: def pseudo_dec(G,L,V,Ord)
365: {
366: N = length(L);
367: S = vector(N);
368: Q = vector(N);
369: R = vector(N);
370: L0 = map(first_element,L);
371: for ( I = 0; I < N; I++ ) {
372: LI = setminus(L0,[L0[I]]);
373: PI = ideal_list_intersection(LI,V,Ord);
374: PI = qsort(PI,comp_tdeg);
375: for ( T = PI; T != []; T = cdr(T) )
376: if ( p_nf(car(T),L0[I],V,Ord) ) break;
377: if ( T == [] ) error("separator : cannot happen");
378: SI = sat_ind(G,car(T),V);
379: QI = SI[0];
380: S[I] = car(T)^SI[1];
381: PV = L[I][1];
382: V0 = setminus(V,PV);
383: #if 0
384: GI = fast_gb(QI,append(V0,PV),0,
385: [[Ord,length(V0)],[Ord,length(PV)]]);
386: #else
387: GI = fast_gb(QI,append(V0,PV),0,
388: [[2,length(V0)],[Ord,length(PV)]]);
389: #endif
390: LCFI = lcfactor(GI,V0,Ord);
391: for ( F = 1, T = LCFI, Gt = QI; T != []; T = cdr(T) ) {
392: St = sat_ind(Gt,T[0],V);
393: Gt = St[0]; F *= T[0]^St[1];
394: }
395: Q[I] = Gt;
396: R[I] = fast_gb(cons(F,QI),V,0,Ord);
397: }
398: return [vtol(Q),vtol(R),vtol(S)];
399: }
400:
401: def iso_comp(G,L,V,Ord)
402: {
403: N = length(L);
404: S = vector(N);
405: Ind = vector(N);
406: Q = vector(N);
407: L0 = map(first_element,L);
1.5 ! noro 408: G = nd_gr(G,V,0,Ord);
1.1 noro 409: for ( I = 0; I < N; I++ ) {
410: LI = setminus(L0,[L0[I]]);
411: PI = ideal_list_intersection(LI,V,Ord);
412: for ( T = PI; T != []; T = cdr(T) )
413: if ( p_nf(car(T),L0[I],V,Ord) ) break;
414: if ( T == [] ) error("separator : cannot happen");
415: S[I] = car(T);
1.5 ! noro 416: QI = sat(G,S[I],V|isgb=1);
1.1 noro 417: PV = L[I][1];
418: V0 = setminus(V,PV);
419: GI = elim_gb(QI,V0,PV,0,[[0,length(V0)],[0,length(PV)]]);
420: Q[I] = [contraction(GI,V0),L0[I]];
421: }
422: return vtol(Q);
423: }
424:
425: /* GTZ primary decompsition */
426:
427: def prima_dec(B,V)
428: {
429: G = nd_gr_trace(B,V,1,GBCheck,0);
430: G0 = G;
431: IntP = [1];
432: QD = [];
433: while ( 1 ) {
434: if ( ideal_inclusion(IntP,G0,V,0) )
435: return QD;
436: W = maxindep(G,V,0); NP = length(W);
437: V0 = setminus(V,W); N = length(V0);
438: V1 = append(V0,W);
439: G1 = fast_gb(G,V1,0,[[0,N],[0,NP]]);
440: LCF = lcfactor(G1,V0,0);
441: L = zprimacomp(G,V0);
442: F = 1;
443: for ( T = LCF, G2 = G1; T != []; T = cdr(T) ) {
444: S = sat_ind(G2,T[0],V1);
445: G2 = S[0]; F *= T[0]^S[1];
446: }
447: for ( T = L, QL = []; T != []; T = cdr(T) )
448: QL = cons(car(T)[0],QL);
449: Int = ideal_list_intersection(QL,V,0);
450: IntP = ideal_intersection(IntP,Int,V,0);
451: QD = append(QD,L);
452: F = p_nf(F,G,V,0);
453: G = cons(F,G);
454: }
455: }
456:
457: /* SL prime decomposition */
458:
459: def prime_dec(B,V)
460: {
461: if ( type(Indep=getopt(indep)) == -1 ) Indep = 0;
462: if ( type(NoLexDec=getopt(nolexdec)) == -1 ) NoLexDec = 0;
463: B = map(sq,B);
464: if ( !NoLexDec )
465: PD = lex_predec1(B,V);
466: else
467: PD = [B];
468: G = ideal_list_intersection(PD,V,0);
469: PD = remove_redundant_comp(G,[],PD,V,0);
470: R = [];
471: for ( T = PD; T != []; T = cdr(T) )
472: R = append(prime_dec_main(car(T),V|indep=Indep),R);
473: if ( Indep ) {
474: G = ideal_list_intersection(map(first_element,R),V,0);
1.5 ! noro 475: if ( !NoLexDec ) R = remove_redundant_comp_first(G,R,V,0);
1.1 noro 476: } else {
477: G = ideal_list_intersection(R,V,0);
1.5 ! noro 478: if ( !NoLexDec ) R = remove_redundant_comp(G,[],R,V,0);
1.1 noro 479: }
480: return R;
481: }
482:
483: def prime_dec_main(B,V)
484: {
485: if ( type(Indep=getopt(indep)) == -1 ) Indep = 0;
486: G = nd_gr_trace(B,V,1,GBCheck,0);
487: IntP = [1];
488: PD = [];
489: while ( 1 ) {
490: /* rad(G) subset IntP */
491: /* check if IntP subset rad(G) */
492: for ( T = IntP; T != []; T = cdr(T) ) {
493: if ( (GNV = modular_radical_membership(car(T),G,V)) ) {
494: F = car(T);
495: break;
496: }
497: }
498: if ( T == [] ) return PD;
499:
500: /* GNV = [GB(<NV*F-1,G>),NV] */
501: G1 = nd_gr_trace(GNV[0],cons(GNV[1],V),1,GBCheck,[[0,1],[0,length(V)]]);
502: G0 = elimination(G1,V);
503: PD0 = zprimecomp(G0,V,Indep);
504: if ( Indep ) {
505: Int = ideal_list_intersection(PD0[0],V,0);
506: IndepSet = PD0[1];
507: for ( PD1 = [], T = PD0[0]; T != []; T = cdr(T) )
508: PD1 = cons([car(T),IndepSet],PD1);
509: PD = append(PD,reverse(PD1));
510: } else {
511: Int = ideal_list_intersection(PD0,V,0);
512: PD = append(PD,PD0);
513: }
514: IntP = ideal_intersection(IntP,Int,V,0);
515: }
516: }
517:
518: /* pre-decomposition */
519:
520: def lex_predec1(B,V)
521: {
522: G = nd_gr_trace(B,V,1,GBCheck,2);
523: for ( T = G; T != []; T = cdr(T) ) {
524: F = fctr(car(T));
525: if ( length(F) > 2 || length(F) == 2 && F[1][1] > 1 ) {
526: for ( R = [], S = cdr(F); S != []; S = cdr(S) ) {
527: Ft = car(S)[0];
528: Gt = map(ptozp,map(p_nf,G,[Ft],V,0));
529: R1 = nd_gr_trace(cons(Ft,Gt),V,1,GBCheck,0);
530: R = cons(R1,R);
531: }
532: return R;
533: }
534: }
535: return [G];
536: }
537:
538: /* zero-dimensional prime/primary decomosition */
539:
540: def zprimedec(B,V,Mod)
541: {
542: L = partial_decomp(B,V,Mod);
543: P = L[0]; NP = L[1];
544: R = [];
545: for ( ; P != []; P = cdr(P) ) R = cons(car(car(P)),R);
546: for ( T = NP; T != []; T = cdr(T) ) {
547: R1 = complete_decomp(car(T),V,Mod);
548: R = append(R1,R);
549: }
550: return R;
551: }
552:
553: def zprimadec(B,V,Mod)
554: {
555: L = partial_qdecomp(B,V,Mod);
556: Q = L[0]; NQ = L[1];
557: R = [];
558: for ( ; Q != []; Q = cdr(Q) ) {
559: T = car(Q); R = cons([T[0],T[1]],R);
560: }
561: for ( T = NQ; T != []; T = cdr(T) ) {
562: R1 = complete_qdecomp(car(T),V,Mod);
563: R = append(R1,R);
564: }
565: return R;
566: }
567:
568: def complete_qdecomp(GD,V,Mod)
569: {
570: GQ = GD[0]; GP = GD[1]; D = GD[2];
571: W = vars(GP);
572: PV = setminus(W,V);
573: N = length(V); PN = length(PV);
574: U = find_npos([GP,D],V,PV,Mod);
575: NV = ttttt;
576: M = gen_minipoly(cons(NV-U,GQ),cons(NV,V),PV,0,NV,Mod);
577: M = ppart(M,NV,Mod);
578: MF = Mod ? modfctr(M) : fctr(M);
579: R = [];
580: for ( T = cdr(MF); T != []; T = cdr(T) ) {
581: S = car(T);
582: Mt = subst(S[0],NV,U);
583: GP1 = fast_gb(cons(Mt,GP),W,Mod,0);
584: GQ1 = fast_gb(cons(Mt^S[1],GQ),W,Mod,0);
585: if ( PV != [] ) {
586: GP1 = elim_gb(GP1,V,PV,Mod,[[0,N],[0,PN]]);
587: GQ1 = elim_gb(GQ1,V,PV,Mod,[[0,N],[0,PN]]);
588: }
589: R = cons([GQ1,GP1],R);
590: }
591: return R;
592: }
593:
594: def partial_qdecomp(B,V,Mod)
595: {
596: Elim = (Elim=getopt(elim))&&type(Elim)!=-1 ? 1 : 0;
597: N = length(V);
598: W = vars(B);
599: PV = setminus(W,V);
600: NP = length(PV);
601: W = append(V,PV);
602: if ( Elim && PV != [] ) Ord = [[0,N],[0,NP]];
603: else Ord = 0;
604: if ( Mod )
605: B = nd_f4(B,W,Mod,Ord);
606: else
607: B = nd_gr_trace(B,W,1,GBCheck,Ord);
608: Q = []; NQ = [[B,B,vector(N+1)]];
609: for ( I = length(V)-1; I >= 0; I-- ) {
610: NQ1 = [];
611: for ( T = NQ; T != []; T = cdr(T) ) {
612: L = partial_qdecomp0(car(T),V,PV,Ord,I,Mod);
613: Q = append(L[0],Q);
614: NQ1 = append(L[1],NQ1);
615: }
616: NQ = NQ1;
617: }
618: return [Q,NQ];
619: }
620:
621: def partial_qdecomp0(GD,V,PV,Ord,I,Mod)
622: {
623: GQ = GD[0]; GP = GD[1]; D = GD[2];
624: N = length(V); PN = length(PV);
625: W = append(V,PV);
626: VI = V[I];
627: M = gen_minipoly(GQ,V,PV,Ord,VI,Mod);
628: M = ppart(M,VI,Mod);
629: if ( Mod )
630: MF = modfctr(M,Mod);
631: else
632: MF = fctr(M);
633: Q = []; NQ = [];
634: if ( length(MF) == 2 && MF[1][1] == 1 ) {
635: D1 = D*1; D1[I] = M;
636: GQelim = elim_gb(GQ,V,PV,Mod,Ord);
637: GPelim = elim_gb(GP,V,PV,Mod,Ord);
638: LD = ldim(GQelim,V);
639: if ( deg(M,VI) == LD )
640: Q = cons([GQelim,GPelim,D1],Q);
641: else
642: NQ = cons([GQelim,GPelim,D1],NQ);
643: return [Q,NQ];
644: }
645: for ( T = cdr(MF); T != []; T = cdr(T) ) {
646: S = car(T); Mt = S[0]; D1 = D*1; D1[I] = Mt;
647:
648: GQ1 = fast_gb(cons(Mt^S[1],GQ),W,Mod,Ord);
649: GQelim = elim_gb(GQ1,V,PV,Mod,Ord);
650: GP1 = fast_gb(cons(Mt,GP),W,Mod,Ord);
651: GPelim = elim_gb(GP1,V,PV,Mod,Ord);
652:
653: D1[N] = LD = ldim(GPelim,V);
654:
655: for ( J = 0; J < N; J++ )
656: if ( D1[J] && deg(D1[J],V[J]) == LD ) break;
657: if ( J < N )
658: Q = cons([GQelim,GPelim,D1],Q);
659: else
660: NQ = cons([GQelim,GPelim,D1],NQ);
661: }
662: return [Q,NQ];
663: }
664:
665: def complete_decomp(GD,V,Mod)
666: {
667: G = GD[0]; D = GD[1];
668: W = vars(G);
669: PV = setminus(W,V);
670: N = length(V); PN = length(PV);
671: U = find_npos(GD,V,PV,Mod);
672: NV = ttttt;
673: M = gen_minipoly(cons(NV-U,G),cons(NV,V),PV,0,NV,Mod);
674: M = ppart(M,NV,Mod);
675: MF = Mod ? modfctr(M) : fctr(M);
676: if ( length(MF) == 2 ) return [G];
677: R = [];
678: for ( T = cdr(MF); T != []; T = cdr(T) ) {
679: Mt = subst(car(car(T)),NV,U);
680: G1 = fast_gb(cons(Mt,G),W,Mod,0);
681: if ( PV != [] ) G1 = elim_gb(G1,V,PV,Mod,[[0,N],[0,PN]]);
682: R = cons(G1,R);
683: }
684: return R;
685: }
686:
687: def partial_decomp(B,V,Mod)
688: {
689: Elim = (Elim=getopt(elim))&&type(Elim)!=-1 ? 1 : 0;
690: N = length(V);
691: W = vars(B);
692: PV = setminus(W,V);
693: NP = length(PV);
694: W = append(V,PV);
695: if ( Elim && PV != [] ) Ord = [[0,N],[0,NP]];
696: else Ord = 0;
697: if ( Mod )
698: B = nd_f4(B,W,Mod,Ord);
699: else
700: B = nd_gr_trace(B,W,1,GBCheck,Ord);
701: P = []; NP = [[B,vector(N+1)]];
702: for ( I = length(V)-1; I >= 0; I-- ) {
703: NP1 = [];
704: for ( T = NP; T != []; T = cdr(T) ) {
705: L = partial_decomp0(car(T),V,PV,Ord,I,Mod);
706: P = append(L[0],P);
707: NP1 = append(L[1],NP1);
708: }
709: NP = NP1;
710: }
711: return [P,NP];
712: }
713:
714: def partial_decomp0(GD,V,PV,Ord,I,Mod)
715: {
716: G = GD[0]; D = GD[1];
717: N = length(V); PN = length(PV);
718: W = append(V,PV);
719: VI = V[I];
720: M = gen_minipoly(G,V,PV,Ord,VI,Mod);
721: M = ppart(M,VI,Mod);
722: if ( Mod )
723: MF = modfctr(M,Mod);
724: else
725: MF = fctr(M);
726: if ( length(MF) == 2 && MF[1][1] == 1 ) {
727: D1 = D*1;
728: D1[I] = M;
729: Gelim = elim_gb(G,V,PV,Mod,Ord);
730: D1[N] = LD = ldim(Gelim,V);
731: GD1 = [Gelim,D1];
732: for ( J = 0; J < N; J++ )
733: if ( D1[J] && deg(D1[J],V[J]) == LD )
734: return [[GD1],[]];
735: return [[],[GD1]];
736: }
737: P = []; NP = [];
738: GI = elim_gb(G,V,PV,Mod,Ord);
739: for ( T = cdr(MF); T != []; T = cdr(T) ) {
740: Mt = car(car(T));
741: D1 = D*1;
742: D1[I] = Mt;
743: GIt = map(p_nf,GI,[Mt],V,Ord);
744: G1 = cons(Mt,GIt);
745: Gelim = elim_gb(G1,V,PV,Mod,Ord);
746: D1[N] = LD = ldim(Gelim,V);
747: for ( J = 0; J < N; J++ )
748: if ( D1[J] && deg(D1[J],V[J]) == LD ) break;
749: if ( J < N )
750: P = cons([Gelim,D1],P);
751: else
752: NP = cons([Gelim,D1],NP);
753: }
754: return [P,NP];
755: }
756:
757: /* prime/primary components over rational function field */
758:
759: def zprimacomp(G,V) {
760: L = zprimadec(G,V,0);
761: R = [];
762: dp_ord(0);
763: for ( T = L; T != []; T = cdr(T) ) {
764: S = car(T);
765: UQ = contraction(S[0],V);
766: UP = contraction(S[1],V);
767: R = cons([UQ,UP],R);
768: }
769: return R;
770: }
771:
772: def zprimecomp(G,V,Indep) {
773: W = maxindep(G,V,0);
774: V0 = setminus(V,W);
775: V1 = append(V0,W);
776: #if 0
777: O1 = [[0,length(V0)],[0,length(W)]];
778: G1 = nd_gr_trace(G,V1,1,GBCheck,O1);
779: dp_ord(0);
780: #else
781: G1 = G;
782: #endif
783: PD = zprimedec(G1,V0,0);
784: dp_ord(0);
785: R = [];
786: for ( T = PD; T != []; T = cdr(T) ) {
787: U = contraction(car(T),V0);
788: R = cons(U,R);
789: }
790: if ( Indep ) return [R,W];
791: else return R;
792: }
793:
794: def fast_gb(B,V,Mod,Ord)
795: {
796: NoRA = (NoRA=getopt(nora))&&type(NoRA)!=-1 ? 1 : 0;
797: if ( Mod )
798: G = nd_f4(B,V,Mod,Ord|nora=NoRA);
799: else {
800: if ( F4 )
801: G = map(ptozp,f4_chrem(B,V,Ord));
802: else
803: G = nd_gr_trace(B,V,1,GBCheck,Ord|nora=NoRA);
804: }
805: return G;
806: }
807:
1.5 ! noro 808: def incremental_gb(A,V,Ord)
! 809: {
! 810: if ( type(Mod=getopt(mod)) == -1 ) Mod = 0;
! 811: if ( type(Block=getopt(gbblock)) == -1 ) Block = 0;
! 812: if ( Mod )
! 813: G = nd_gr(A,V,Mod,Ord);
! 814: else if ( Procs ) {
! 815: Arg0 = ["nd_gr",A,V,0,Ord];
! 816: Arg1 = ["nd_gr_trace",A,V,1,GBCheck,Ord];
! 817: G = competitive_exec(Procs,Arg0,Arg1);
! 818: } else if ( Block )
! 819: G = nd_gr(A,V,0,Ord|gbblock=Block);
! 820: else
! 821: G = nd_gr(A,V,0,Ord);
! 822: return G;
! 823: }
1.1 noro 824:
825: def elim_gb(G,V,PV,Mod,Ord)
826: {
827: N = length(V); PN = length(PV);
828: O1 = [[0,N],[0,PN]];
829: if ( Ord == O1 )
830: Ord = Ord[0][0];
831: if ( Mod ) /* XXX */
832: G = dp_gr_mod_main(G,V,0,Mod,Ord);
833: else if ( Procs ) {
834: Arg0 = ["nd_gr_trace",G,V,1,GBCheck,Ord];
835: Arg1 = ["nd_gr_trace_rat",G,V,PV,1,GBCheck,O1,Ord];
836: G = competitive_exec(Procs,Arg0,Arg1);
837: } else
838: G = nd_gr_trace(G,V,1,GBCheck,Ord);
839: return G;
840: }
841:
842: def ldim(G,V)
843: {
844: O0 = dp_ord(); dp_ord(0);
845: D = length(dp_mbase(map(dp_ptod,G,V)));
846: dp_ord(O0);
847: return D;
848: }
849:
850: def make_mod_subst(GD,V,PV,HC)
851: {
852: N = length(V);
853: PN = length(PV);
854: G = GD[0]; D = GD[1];
855: for ( I = 0; ; I = (I+1)%100 ) {
856: Mod = lprime(I);
857: S = [];
858: for ( J = PN-1; J >= 0; J-- )
859: S = append([PV[J],random()%Mod],S);
860: for ( T = HC; T != []; T = cdr(T) )
861: if ( !(subst(car(T),S)%Mod) ) break;
862: if ( T != [] ) continue;
863: for ( J = 0; J < N; J++ ) {
864: M = subst(D[J],S);
865: F = modsqfr(M,Mod);
866: if ( length(F) != 2 || F[1][1] != 1 ) break;
867: }
868: if ( J < N ) continue;
869: G0 = map(subst,G,S);
870: return [G0,Mod];
871: }
872: }
873:
874: def rsgn()
875: {
876: return random()%2 ? 1 : -1;
877: }
878:
879: def find_npos(GD,V,PV,Mod)
880: {
881: N = length(V); PN = length(PV);
882: G = GD[0]; D = GD[1]; LD = D[N];
883: Ord0 = dp_ord(); dp_ord(0);
884: HC = map(dp_hc,map(dp_ptod,G,V));
885: dp_ord(Ord0);
886: if ( !Mod ) {
887: W = append(V,PV);
888: G1 = nd_gr_trace(G,W,1,GBCheck,[[0,N],[0,PN]]);
889: L = make_mod_subst([G1,D],V,PV,HC);
890: return find_npos([L[0],D],V,[],L[1]);
891: }
892: N = length(V);
893: NV = ttttt;
894: for ( B = 2; ; B++ ) {
895: for ( J = N-2; J >= 0; J-- ) {
896: for ( U = 0, K = J; K < N; K++ )
897: U += rsgn()*((random()%B+1))*V[K];
898: M = minipolym(G,V,0,U,NV,Mod);
899: if ( deg(M,NV) == LD ) return U;
900: }
901: }
902: }
903:
904: def gen_minipoly(G,V,PV,Ord,VI,Mod)
905: {
906: if ( PV == [] ) {
907: NV = ttttt;
908: if ( Mod )
909: M = minipolym(G,V,Ord,VI,NV,Mod);
910: else
911: M = minipoly(G,V,Ord,VI,NV);
912: return subst(M,NV,VI);
913: }
914: W = setminus(V,[VI]);
915: PV1 = cons(VI,PV);
916: #if 0
917: while ( 1 ) {
918: V1 = append(W,PV1);
919: if ( Mod )
920: G = nd_gr(G,V1,Mod,[[0,1],[0,length(V1)-1]]|nora=1);
921: else
922: G = nd_gr_trace(G,V1,1,GBCheck,[[0,1],[0,length(V1)-1]]|nora=1);
923: if ( W == [] ) break;
924: else {
925: W = cdr(W);
926: G = elimination(G,cdr(V1));
927: }
928: }
929: #elif 1
930: if ( Mod ) {
931: G = nd_gr(G,V1,Mod,[[0,length(W)],[0,length(PV1)]]|nora=1);
932: G = elimination(G,PV1);
933: } else {
934: PV2 = setminus(PV1,[PV1[length(PV1)-1]]);
935: V2 = append(W,PV2);
936: G = nd_gr_trace(G,V2,1,GBCheck,[[0,length(W)],[0,length(PV2)]]|nora=1);
937: G = elimination(G,PV1);
938: }
939: #else
940: V1 = append(W,PV1);
941: if ( Mod )
942: G = nd_gr(G,V1,Mod,[[0,length(W)],[0,length(PV1)]]|nora=1);
943: else
944: G = nd_gr_trace(G,V1,1,GBCheck,[[0,length(W)],[0,length(PV1)]]|nora=1);
945: G = elimination(G,PV1);
946: #endif
947: if ( Mod )
948: G = nd_gr(G,PV1,Mod,[[0,1],[0,length(PV)]]|nora=1);
949: else
950: G = nd_gr_trace(G,PV1,1,GBCheck,[[0,1],[0,length(PV)]]|nora=1);
951: for ( M = car(G), T = cdr(G); T != []; T = cdr(T) )
952: if ( deg(car(T),VI) < deg(M,VI) ) M = car(T);
953: return M;
954: }
955:
956: def indepset(V,H)
957: {
958: if ( H == [] ) return V;
959: N = -1;
960: for ( T = V; T != []; T = cdr(T) ) {
961: VI = car(T);
962: HI = [];
963: for ( S = H; S != []; S = cdr(S) )
964: if ( !tdiv(car(S),VI) ) HI = cons(car(S),HI);
965: RI = indepset(setminus(V,[VI]),HI);
966: if ( length(RI) > N ) {
967: R = RI; N = length(RI);
968: }
969: }
970: return R;
971: }
972:
973: def maxindep(B,V,O)
974: {
975: G = nd_gr_trace(B,V,1,GBCheck,O);
976: Old = dp_ord();
977: dp_ord(O);
978: H = map(dp_dtop,map(dp_ht,map(dp_ptod,G,V)),V);
979: H = dp_mono_raddec(H,V);
980: N = length(V);
981: Dep = [];
982: for ( T = H, Len = N+1; T != []; T = cdr(T) ) {
983: M = length(car(T));
984: if ( M < Len ) {
985: Dep = [car(T)];
986: Len = M;
987: } else if ( M == Len )
988: Dep = cons(car(T),Dep);
989: }
990: R = setminus(V,Dep[0]);
991: dp_ord(Old);
992: return R;
993: }
994:
995: /* ideal operations */
996: def contraction(G,V)
997: {
998: C = [];
999: for ( T = G; T != []; T = cdr(T) ) {
1000: C1 = dp_hc(dp_ptod(car(T),V));
1001: S = fctr(C1);
1002: for ( S = cdr(S); S != []; S = cdr(S) )
1003: if ( !member(S[0][0],C) ) C = cons(S[0][0],C);
1004: }
1005: W = vars(G);
1006: PV = setminus(W,V);
1007: W = append(V,PV);
1008: NV = ttttt;
1009: for ( T = C, S = 1; T != []; T = cdr(T) )
1010: S *= car(T);
1011: G = saturation([G,NV],S,W);
1012: return G;
1013: }
1014:
1015: def ideal_list_intersection(L,V,Ord)
1016: {
1017: if ( type(Mod=getopt(mod)) == -1 ) Mod = 0;
1018: N = length(L);
1019: if ( N == 0 ) return [1];
1020: if ( N == 1 ) return fast_gb(L[0],V,Mod,Ord);
1021: N2 = idiv(N,2);
1022: for ( L1 = [], I = 0; I < N2; I++ ) L1 = cons(L[I],L1);
1023: for ( L2 = []; I < N; I++ ) L2 = cons(L[I],L2);
1024: I1 = ideal_list_intersection(L1,V,Ord|mod=Mod);
1025: I2 = ideal_list_intersection(L2,V,Ord|mod=Mod);
1026: return ideal_intersection(I1,I2,V,Ord|mod=Mod,
1027: gbblock=[[0,length(I1)],[length(I1),length(I2)]]);
1028: }
1029:
1030: def ideal_intersection(A,B,V,Ord)
1031: {
1032: if ( type(Mod=getopt(mod)) == -1 ) Mod = 0;
1033: if ( type(Block=getopt(gbblock)) == -1 ) Block = 0;
1034: T = ttttt;
1035: if ( Mod )
1036: G = nd_gr(append(vtol(ltov(A)*T),vtol(ltov(B)*(1-T))),
1037: cons(T,V),Mod,[[0,1],[Ord,length(V)]]);
1038: else
1039: if ( Procs ) {
1040: Arg0 = ["nd_gr",
1041: append(vtol(ltov(A)*T),vtol(ltov(B)*(1-T))),
1042: cons(T,V),0,[[0,1],[Ord,length(V)]]];
1043: Arg1 = ["nd_gr_trace",
1044: append(vtol(ltov(A)*T),vtol(ltov(B)*(1-T))),
1045: cons(T,V),1,GBCheck,[[0,1],[Ord,length(V)]]];
1046: G = competitive_exec(Procs,Arg0,Arg1);
1047: } else {
1048: if ( Block )
1049: G = nd_gr(append(vtol(ltov(A)*T),vtol(ltov(B)*(1-T))),
1050: cons(T,V),0,[[0,1],[Ord,length(V)]]|gbblock=Block);
1051: else
1052: G = nd_gr(append(vtol(ltov(A)*T),vtol(ltov(B)*(1-T))),
1053: cons(T,V),0,[[0,1],[Ord,length(V)]]);
1054: }
1055: G0 = elimination(G,V);
1056: return G0;
1057: }
1058:
1059: /* returns GB if F notin rad(G) */
1060:
1061: def radical_membership(F,G,V) {
1062: F = p_nf(F,G,V,0);
1063: if ( !F ) return 0;
1064: NV = ttttt;
1065: T = nd_gr_trace(cons(NV*F-1,G),cons(NV,V),1,GBCheck,0);
1066: if ( type(car(T)) != 1 ) return [T,NV];
1067: else return 0;
1068: }
1069:
1070: def quick_radical_membership(F,G,V) {
1071: F = p_nf(F,G,V,0);
1072: if ( !F ) return 1;
1073: NV = ttttt;
1074: T = nd_f4(cons(NV*F-1,G),cons(NV,V),lprime(0),0);
1075: if ( type(car(T)) != 1 ) return 0;
1076: else return 1;
1077: }
1078:
1079: def modular_radical_membership(F,G,V) {
1080: F = p_nf(F,G,V,0);
1081: if ( !F ) return 0;
1082: NV = ttttt;
1083: for ( J = 0; ; J++ ) {
1084: Mod = lprime(J);
1085: H = map(dp_hc,map(dp_ptod,G,V));
1086: for ( ; H != []; H = cdr(H) ) if ( !(car(H)%Mod) ) break;
1087: if ( H != [] ) continue;
1088:
1089: T = nd_f4(cons(NV*F-1,G),cons(NV,V),Mod,0);
1090: if ( type(car(T)) == 1 ) {
1091: I = radical_membership_rep(F,G,V,-1,0,Mod);
1092: I1 = radical_membership_rep(F,G,V,I,0,0);
1093: if ( I1 > 0 ) return 0;
1094: }
1095: return radical_membership(F,G,V);
1096: }
1097: }
1098:
1099: def radical_membership_rep(F,G,V,Max,Ord,Mod) {
1100: Ft = F;
1101: I = 1;
1102: while ( Max < 0 || I <= Max ) {
1103: Ft = nd_nf(Ft,G,V,Ord,Mod);
1104: if ( !Ft ) return I;
1105: Ft *= F;
1106: I++;
1107: }
1108: return -1;
1109: }
1110:
1111: def ideal_product(A,B,V)
1112: {
1113: dp_ord(0);
1114: DA = map(dp_ptod,A,V);
1115: DB = map(dp_ptod,B,V);
1116: DegA = map(dp_td,DA);
1117: DegB = map(dp_td,DB);
1118: for ( PA = [], T = A, DT = DegA; T != []; T = cdr(T), DT = cdr(DT) )
1119: PA = cons([car(T),car(DT)],PA);
1120: PA = reverse(PA);
1121: for ( PB = [], T = B, DT = DegB; T != []; T = cdr(T), DT = cdr(DT) )
1122: PB = cons([car(T),car(DT)],PB);
1123: PB = reverse(PB);
1124: R = [];
1125: for ( T = PA; T != []; T = cdr(T) )
1126: for ( S = PB; S != []; S = cdr(S) )
1127: R = cons([car(T)[0]*car(S)[0],car(T)[1]+car(S)[1]],R);
1128: T = qsort(R,comp_by_second);
1129: T = map(first_element,T);
1130: Len = length(A)>length(B)?length(A):length(B);
1131: Len *= 2;
1132: L = sep_list(T,Len); B0 = L[0]; B1 = L[1];
1133: R = nd_gr_trace(B0,V,0,-1,0);
1134: while ( B1 != [] ) {
1135: print(length(B1));
1136: L = sep_list(B1,Len);
1137: B0 = L[0]; B1 = L[1];
1138: R = nd_gr_trace(append(R,B0),V,0,-1,0|gbblock=[[0,length(R)]],nora=1);
1139: }
1140: return R;
1141: }
1142:
1143: def saturation(GNV,F,V)
1144: {
1145: G = GNV[0]; NV = GNV[1];
1146: if ( Procs ) {
1147: Arg0 = ["nd_gr_trace",
1148: cons(NV*F-1,G),cons(NV,V),0,GBCheck,[[0,1],[0,length(V)]]];
1149: Arg1 = ["nd_gr_trace",
1150: cons(NV*F-1,G),cons(NV,V),1,GBCheck,[[0,1],[0,length(V)]]];
1151: G1 = competitive_exec(Procs,Arg0,Arg1);
1152: } else
1153: G1 = nd_gr_trace(cons(NV*F-1,G),cons(NV,V),SatHomo,GBCheck,[[0,1],[0,length(V)]]);
1154: return elimination(G1,V);
1155: }
1156:
1157: def sat(G,F,V)
1158: {
1.5 ! noro 1159: if ( type(IsGB=getopt(isgb)) == -1 ) IsGB = 0;
1.1 noro 1160: NV = ttttt;
1161: if ( Procs ) {
1162: Arg0 = ["nd_gr_trace",
1163: cons(NV*F-1,G),cons(NV,V),0,GBCheck,[[0,1],[0,length(V)]]];
1164: Arg1 = ["nd_gr_trace",
1165: cons(NV*F-1,G),cons(NV,V),1,GBCheck,[[0,1],[0,length(V)]]];
1166: G1 = competitive_exec(Procs,Arg0,Arg1);
1.5 ! noro 1167: } else {
! 1168: B1 = append(G,[NV*F-1]);
! 1169: V1 = cons(NV,V);
! 1170: Ord1 = [[0,1],[0,length(V)]];
! 1171: if ( IsGB )
! 1172: G1 = nd_gr_trace(B1,V1,SatHomo,GBCheck,Ord1|
! 1173: gbblock=[[0,length(G)]]);
! 1174: else
! 1175: G1 = nd_gr_trace(B1,V1,SatHomo,GBCheck,Ord1);
! 1176: }
1.1 noro 1177: return elimination(G1,V);
1178: }
1179:
1180: def satind(G,F,V)
1181: {
1182: NV = ttttt;
1183: N = length(V);
1184: B = append(G,[NV*F-1]);
1185: V1 = cons(NV,V);
1186: D = nd_gr_trace(B,V1,1,GBCheck,[[0,1],[0,N]]
1187: |nora=1,gentrace=1,gbblock=[[0,length(G)]]);
1188: G1 = D[0];
1189: Len = length(G1);
1190: Deg = compute_deg(B,V1,NV,D);
1191: D1 = 0;
1192: R = [];
1193: M = length(B);
1194: for ( I = 0; I < Len; I++ ) {
1195: if ( !member(NV,vars(G1[I])) ) {
1196: for ( J = 1; J < M; J++ )
1197: D1 = MAX(D1,Deg[I][J]);
1198: R = cons(G1[I],R);
1199: }
1200: }
1201: return [reverse(R),D1];
1202: }
1203:
1204: def sat_ind(G,F,V)
1205: {
1206: NV = ttttt;
1207: F = p_nf(F,G,V,0);
1208: for ( I = 0, GI = G; ; I++ ) {
1209: G1 = colon(GI,F,V);
1210: if ( ideal_inclusion(G1,GI,V,0) ) {
1211: return [GI,I];
1212: }
1213: else GI = G1;
1214: }
1215: }
1216:
1217: def colon(G,F,V)
1218: {
1.5 ! noro 1219: if ( type(IsGB=getopt(isgb)) == -1 ) IsGB = 0;
1.1 noro 1220: F = p_nf(F,G,V,0);
1221: if ( !F ) return [1];
1.5 ! noro 1222: if ( IsGB )
! 1223: T = ideal_intersection(G,[F],V,0|gbblock=[[0,length(G)]]);
! 1224: else
! 1225: T = ideal_intersection(G,[F],V,0);
1.1 noro 1226: return map(ptozp,map(sdiv,T,F));
1227: }
1228:
1229: def ideal_colon(G,F,V)
1230: {
1231: G = nd_gr(G,V,0,0);
1.5 ! noro 1232: for ( T = F, L = []; T != []; T = cdr(T) )
! 1233: L = cons(colon(G,car(T),V|isgb=1),L);
! 1234: L = reverse(L);
1.1 noro 1235: return ideal_list_intersection(L,V,0);
1236: }
1237:
1.2 noro 1238: def ideal_sat(G,F,V)
1239: {
1240: G = nd_gr(G,V,0,0);
1241: L = mapat(sat,1,G,F,V);
1242: return ideal_list_intersection(L,V,0);
1243: }
1244:
1.1 noro 1245: def ideal_inclusion(F,G,V,O)
1246: {
1247: if ( type(Mod=getopt(mod)) == -1 ) Mod = 0;
1248: if ( Mod ) {
1249: for ( T = F; T != []; T = cdr(T) )
1250: if ( p_nf_mod(car(T),G,V,O,Mod) ) return 0;
1251: } else {
1252: for ( T = F; T != []; T = cdr(T) )
1253: if ( p_nf(car(T),G,V,O) ) return 0;
1254: }
1255: return 1;
1256: }
1257:
1258: /* remove redundant components */
1259:
1260: def qd_simp_comp(QP,V)
1261: {
1262: R = ltov(QP);
1263: N = length(R);
1264: for ( I = 0; I < N; I++ ) {
1265: if ( R[I] ) {
1266: QI = R[I][0]; PI = R[I][1];
1267: for ( J = I+1; J < N; J++ )
1268: if ( R[J] && gb_comp(PI,R[J][1]) ) {
1269: QI = ideal_intersection(QI,R[J][0],V,0);
1270: R[J] = 0;
1271: }
1272: R[I] = [QI,PI];
1273: }
1274: }
1275: for ( I = N-1, S = []; I >= 0; I-- )
1276: if ( R[I] ) S = cons(R[I],S);
1277: return S;
1278: }
1279:
1280: def qd_remove_redundant_comp(G,Iso,Emb,V,Ord)
1281: {
1282: IsoInt = ideal_list_intersection(map(first_element,Iso),V,Ord);
1283: Emb = qd_simp_comp(Emb,V);
1284: Emb = qsort(Emb);
1285: A = ltov(Emb);
1286: N = length(A);
1287: for ( I = 0; I < N; I++ ) {
1288: if ( !A[I] ) continue;
1289: for ( T = [], J = 0; J < N; J++ )
1290: if ( J != I && A[J] ) T = cons(A[J][0],T);
1291: Int = ideal_list_intersection(T,V,Ord);
1292: Int = ideal_intersection(IsoInt,Int,V,Ord);
1293: if ( gb_comp(Int,G) ) A[I] = 0;
1294: }
1295: for ( T = [], I = 0; I < N; I++ )
1296: if ( A[I] ) T = cons(A[I],T);
1297: return reverse(T);
1298: }
1299:
1300: def remove_redundant_comp(G,Iso,Emb,V,Ord)
1301: {
1302: IsoInt = ideal_list_intersection(Iso,V,Ord);
1303:
1304: A = ltov(Emb);
1305: N = length(A);
1306: for ( I = 0; I < N; I++ ) {
1307: if ( !A[I] ) continue;
1308: for ( J = I+1; J < N; J++ )
1309: if ( A[J] && gb_comp(A[I],A[J]) ) A[J] = 0;
1310: }
1311: for ( I = 0; I < N; I++ ) {
1312: if ( !A[I] ) continue;
1313: for ( T = [], J = 0; J < N; J++ )
1314: if ( J != I && A[J] ) T = cons(A[J],T);
1315: Int = ideal_list_intersection(cons(IsoInt,T),V,Ord);
1316: if ( gb_comp(Int,G) ) A[I] = 0;
1317: }
1318: for ( T = [], I = 0; I < N; I++ )
1319: if ( A[I] ) T = cons(A[I],T);
1320: return reverse(T);
1321: }
1322:
1323: def remove_redundant_comp_first(G,P,V,Ord)
1324: {
1325: A = ltov(P);
1326: N = length(A);
1327: for ( I = 0; I < N; I++ ) {
1328: if ( !A[I] ) continue;
1329: for ( J = I+1; J < N; J++ )
1330: if ( A[J] && gb_comp(A[I][0],A[J][0]) ) A[J] = 0;
1331: }
1332: for ( I = 0; I < N; I++ ) {
1333: if ( !A[I] ) continue;
1334: for ( T = [], J = 0; J < N; J++ )
1335: if ( J != I && A[J] ) T = cons(A[J][0],T);
1336: Int = ideal_list_intersection(T,V,Ord);
1337: if ( gb_comp(Int,G) ) A[I] = 0;
1338: }
1339: for ( T = [], I = 0; I < N; I++ )
1340: if ( A[I] ) T = cons(A[I],T);
1341: return reverse(T);
1342: }
1343:
1344: /* polynomial operations */
1345:
1346: def ppart(F,V,Mod)
1347: {
1348: if ( !Mod )
1349: G = nd_gr([F],[V],0,0);
1350: else
1351: G = dp_gr_mod_main([F],[V],0,Mod,0);
1352: return G[0];
1353: }
1354:
1355:
1356: def sq(F)
1357: {
1358: if ( !F ) return 0;
1359: A = cdr(fctr(F));
1360: for ( R = 1; A != []; A = cdr(A) )
1361: R *= car(car(A));
1362: return R;
1363: }
1364:
1365: def lcfactor(G,V,O)
1366: {
1367: O0 = dp_ord(); dp_ord(O);
1368: C = [];
1369: for ( T = G; T != []; T = cdr(T) ) {
1370: C1 = dp_hc(dp_ptod(car(T),V));
1371: S = fctr(C1);
1372: for ( S = cdr(S); S != []; S = cdr(S) )
1373: if ( !member(S[0][0],C) ) C = cons(S[0][0],C);
1374: }
1375: dp_ord(O0);
1376: return C;
1377: }
1378:
1379: /* Ti = [D,I,M,C] */
1380:
1381: def compute_deg0(Ti,P,V,TV)
1382: {
1383: N = length(P[0]);
1384: Num = vector(N);
1385: for ( I = 0; I < N; I++ ) Num[I] = -1;
1386: for ( ; Ti != []; Ti = cdr(Ti) ) {
1387: Sj = car(Ti);
1388: Dj = Sj[0];
1389: Ij =Sj[1];
1390: Mj = deg(type(Sj[2])==9?dp_dtop(Sj[2],V):Sj[2],TV);
1391: Pj = P[Ij];
1392: if ( Dj )
1393: for ( I = 0; I < N; I++ )
1394: if ( Pj[I] >= 0 ) {
1395: T = Mj+Pj[I];
1396: Num[I] = MAX(Num[I],T);
1397: }
1398: }
1399: return Num;
1400: }
1401:
1402: def compute_deg(B,V,TV,Data)
1403: {
1404: GB = Data[0];
1405: Homo = Data[1];
1406: Trace = Data[2];
1407: IntRed = Data[3];
1408: Ind = Data[4];
1409: DB = map(dp_ptod,B,V);
1410: if ( Homo ) {
1411: DB = map(dp_homo,DB);
1412: V0 = append(V,[hhh]);
1413: } else
1414: V0 = V;
1415: Perm = Trace[0]; Trace = cdr(Trace);
1416: for ( I = length(Perm)-1, T = Trace; T != []; T = cdr(T) )
1417: if ( (J=car(T)[0]) > I ) I = J;
1418: N = I+1;
1419: N0 = length(B);
1420: P = vector(N);
1421: for ( T = Perm, I = 0; T != []; T = cdr(T), I++ ) {
1422: Pi = car(T);
1423: C = vector(N0);
1424: for ( J = 0; J < N0; J++ ) C[J] = -1;
1425: C[Pi[1]] = 0;
1426: P[Pi[0]] = C;
1427: }
1428: for ( T = Trace; T != []; T = cdr(T) ) {
1429: Ti = car(T); P[Ti[0]] = compute_deg0(Ti[1],P,V0,TV);
1430: }
1431: M = length(Ind);
1432: for ( T = IntRed; T != []; T = cdr(T) ) {
1433: Ti = car(T); P[Ti[0]] = compute_deg0(Ti[1],P,V,TV);
1434: }
1435: R = [];
1436: for ( J = 0; J < M; J++ ) {
1437: U = P[Ind[J]];
1438: R = cons(U,R);
1439: }
1440: return reverse(R);
1441: }
1442:
1443: /* set theoretic functions */
1444:
1445: def member(A,S)
1446: {
1447: for ( ; S != []; S = cdr(S) )
1448: if ( car(S) == A ) return 1;
1449: return 0;
1450: }
1451:
1452: def elimination(G,V) {
1453: for ( R = [], T = G; T != []; T = cdr(T) )
1454: if ( setminus(vars(car(T)),V) == [] ) R =cons(car(T),R);
1455: return R;
1456: }
1457:
1458: def setintersection(A,B)
1459: {
1460: for ( L = []; A != []; A = cdr(A) )
1461: if ( member(car(A),B) )
1462: L = cons(car(A),L);
1463: return L;
1464: }
1465:
1466: def setminus(A,B) {
1467: for ( T = reverse(A), R = []; T != []; T = cdr(T) ) {
1468: for ( S = B, M = car(T); S != []; S = cdr(S) )
1469: if ( car(S) == M ) break;
1470: if ( S == [] ) R = cons(M,R);
1471: }
1472: return R;
1473: }
1474:
1475: def sep_list(L,N)
1476: {
1477: if ( length(L) <= N ) return [L,[]];
1478: R = [];
1479: for ( T = L, I = 0; I < N; I++, T = cdr(T) )
1480: R = cons(car(T),R);
1481: return [reverse(R),T];
1482: }
1483:
1484: def first_element(L)
1485: {
1486: return L[0];
1487: }
1488:
1489: def comp_tdeg(A,B)
1490: {
1491: DA = tdeg(A);
1492: DB = tdeg(B);
1493: if ( DA > DB ) return 1;
1494: else if ( DA < DB ) return -1;
1495: else return 0;
1496: }
1497:
1498: def tdeg(P)
1499: {
1500: dp_ord(0);
1501: return dp_td(dp_ptod(P,vars(P)));
1502: }
1503:
1504: def comp_by_ord(A,B)
1505: {
1506: if ( dp_ht(A) > dp_ht(B) ) return 1;
1507: else if ( dp_ht(A) < dp_ht(B) ) return -1;
1508: else return 0;
1509: }
1510:
1511: def comp_by_second(A,B)
1512: {
1513: if ( A[1] > B[1] ) return 1;
1514: else if ( A[1] < B[1] ) return -1;
1515: else return 0;
1516: }
1517: endmodule$
1518: end$