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