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