Annotation of OpenXM/src/asir-contrib/testing/noro/module_syz.rr, Revision 1.6
1.1 noro 1: module newsyz;
2:
1.4 noro 3: localf module_syz, module_syz_old;
1.3 noro 4: localf simplify_syz, icont, mod, remove_cont,ordcheck;
5: localf complsb, complsb_sd, sortlsb, find_pos, find_pos, reduce, lres_setup, dpm_sort1, comp_pos;
6: localf fres,minres,sres,minsres,lres, create_base_ord, simplify_k, simplify_by_k, remove_k, remove_k1, extract_nonzero;
7: localf nonzero, phi, syz_check, renumber_pos, compress, compress_h;
8: localf syz_check0,phi0,todpmlist,dpmlisttollist;
1.1 noro 9:
10: /* F : a list of (lists or polynomials),
11: V : a variable list, H >1=> over GF(H), H=0,1=> over Q
12: O : term order
13: return: [GS,G]
14: GS : a GB of syz(F) wrt [1,O] (POT), G: a GB of F wrt [1,O]
15: */
16:
1.3 noro 17: // return [dpmlist F,rank N]
18: def todpmlist(F,V)
1.1 noro 19: {
1.3 noro 20: K = length(F);
21: for ( I = 0; I < K; I++ ) if ( F[I] ) break;
22: if ( I == K ) return [];
23: if ( type(F[I]) <= 2 ) {
24: // F is a polynimial list
25: F = map(dp_ptod,F,V);
26: F = map(dpm_dptodpm,F,1);
27: N = 1;
28: } else if ( type(F[I]) == 9 ) {
29: // F is a dpoly list
30: F = map(dpm_dptodpm,F,1);
31: N = 1;
32: } else if ( type(F[I]) == 4 ) {
33: // F is a list of poly lists
34: N = length(F[0]);
35: F = map(dpm_ltod,F,V);
36: } else if ( type(F[I]) == 26 ) {
37: // F is a DPM list
38: for ( N = 0, T = F; T != []; T = cdr(T) ) {
39: for ( A = car(T); A; A = dpm_rest(A) ) {
40: N1 = dpm_hp(A);
41: if ( N1 > N ) N = N1;
42: }
43: }
44: } else {
45: error("todpmlist: arugument type is invalid.");
46: }
47: return [F,N];
48: }
49:
50: def module_syz(F,V,H,Ord)
51: {
52: if ( type(Weyl=getopt(weyl)) == -1 ) Weyl = 0;
53: if ( type(DP=getopt(dp)) == -1 ) DP = 0;
54: if ( type(F4=getopt(f4)) == -1 ) F4 = 0;
55: dp_ord(Ord);
56: K = length(F);
57: for ( I = 0; I < K; I++ ) if ( F[I] ) break;
58: if ( I == K ) return [[],[],[]];
59: L = todpmlist(F,V);
60: F = L[0]; N = L[1];
61: dp_ord([1,Ord]);
62: B = [];
63: for ( I = 0; I < K; I++ ) {
64: B = cons(F[I]+dpm_dptodpm(dp_ptod(1,V),N+I+1),B);
65: }
66: B = reverse(B);
1.4 noro 67: if ( H >= 2 ) {
68: // finite field
69: if ( Weyl )
70: G = nd_weyl_gr(B,V,H,[1,Ord]|dp=1);
71: else if ( F4 )
72: G = nd_f4(B,V,H,[1,Ord]|dp=1);
73: else
74: G = nd_gr(B,V,H,[1,Ord]|dp=1);
75: } else {
76: if ( Weyl )
77: G = nd_weyl_gr(B,V,0,[1,Ord]|dp=1,homo=H);
78: else if ( F4 ) {
79: Ind = 0;
80: while ( 1 ) {
81: G = nd_f4_trace(B,V,H,-lprime(Ind),[1,Ord]|dp=1);
82: if ( G ) break;
83: else Ind++;
84: }
85: } else
1.6 ! noro 86: G = nd_gr_trace(B,V,H,1,[1,Ord]|dp=1);
1.4 noro 87: }
1.3 noro 88: G0 = []; S0 = []; Gen0 = [];
89: for ( T = G; T != []; T = cdr(T) ) {
90: H = car(T);
91: if ( dpm_hp(H) > N ) {
92: S0 = cons(dpm_shift(H,N),S0);
93: } else {
94: L = dpm_split(H,N);
95: G0 = cons(L[0],G0);
96: Gen0 = cons(dpm_shift(L[1],N),Gen0);
1.1 noro 97: }
1.3 noro 98: }
99: #if 0
100: S0 = reverse(S0); G0 = reverse(G0); Gen0 = reverse(Gen0);
101: #endif
102: if ( !DP ) {
103: S0 = map(dpm_dtol,S0,V); G0 = map(dpm_dtol,G0,V); Gen0 = map(dpm_dtol,Gen0,V);
104: }
105: return [S0,G0,Gen0];
106: }
107:
1.4 noro 108: def module_syz_old(F,V,H,O)
109: {
110: Weyl = type(getopt(weyl)) != -1 ? 1 : 0;
111: K = length(F);
112: if ( type(F[0]) <= 2 ) {
113: for ( T = [], S = F; S != []; S = cdr(S) )
114: T = cons([car(S)],T);
115: F = reverse(T);
116: }
117: N = length(F[0]);
118: B = [];
119: for ( I = 0; I < K; I++ ) {
120: E = vector(N+K);
121: for ( J = 0; J < N; J++ ) E[J] = F[I][J];
122: E[N+I] = 1;
123: B = cons(vtol(E),B);
124: }
125: B = reverse(B);
126: if ( H >= 2 ) {
127: if ( Weyl )
128: G = nd_weyl_gr(B,V,H,[1,O]);
129: else
130: G = nd_gr(B,V,H,[1,O]);
131: } else {
132: if ( Weyl )
133: G = nd_weyl_gr_trace(B,V,H,-1,[1,O]);
134: else
135: G = nd_gr_trace(B,V,H,-1,[1,O]);
136: }
137: G0 = []; S0 = []; Gen0 = [];
138: for ( T = G; T != []; T = cdr(T) ) {
139: H = car(T);
140: for ( J = 0; J < N; J++ ) if ( H[J] ) break;
141: if ( J == N ) {
142: H1 = vector(K);
143: for ( J = 0; J < K; J++ ) H1[J] = H[N+J];
144: S0 = cons(vtol(H1),S0);
145: } else {
146: H1 = vector(N);
147: for ( J = 0; J < N; J++ ) H1[J] = H[J];
148: G0 = cons(vtol(H1),G0);
149: H1 = vector(K);
150: for ( J = 0; J < K; J++ ) H1[J] = H[N+J];
151: Gen0 = cons(vtol(H1),Gen0);
152: }
153: }
154: return [S0,G0,Gen0];
155: }
156:
1.3 noro 157: def fres(F,V,H,O)
158: {
159: if ( type(Weyl=getopt(weyl)) == -1 ) Weyl = 0;
160: if ( type(F4=getopt(f4)) == -1 ) F4 = 0;
161: if ( type(DP=getopt(dp)) == -1 ) DP = 0;
162: L = todpmlist(F,V); F = L[0];
163: R = [F];
164: while ( 1 ) {
165: if ( Weyl )
166: L = module_syz(car(R),V,H,O|weyl=1,dp=1,f4=F4);
167: else
168: L = module_syz(car(R),V,H,O|dp=1,f4=F4);
169: L = L[0];
170: L = ordcheck(L,V);
171: if ( L == [] ) {
1.5 noro 172: R = reverse(R);
1.3 noro 173: if ( DP ) return R;
174: else return map(dpmlisttollist,R,V);
175: }
176: R = cons(L,R);
177: }
178: }
179:
180: def minres(F,V,H,O)
181: {
182: if ( type(Weyl=getopt(weyl)) == -1 ) Weyl = 0;
183: if ( type(F4=getopt(f4)) == -1 ) F4 = 0;
184: if ( type(DP=getopt(dp)) == -1 ) DP = 0;
185: L = todpmlist(F,V); F = L[0];
186: R = [F];
187: while ( 1 ) {
188: if ( Weyl )
189: L = module_syz(car(R),V,H,O|weyl=1,dp=1,f4=F4);
190: else
191: L = module_syz(car(R),V,H,O|dp=1,f4=F4);
192: L = L[0];
193: L = ordcheck(L,V);
194: if ( L == [] ) break;
195: S = dpm_simplify_syz(L,R[0]);
196: if ( S[0] == [] && S[1] == [] ) {
197: R = cdr(R); break;
198: }
199: R = append([S[0],S[1]],cdr(R));
200: if ( R[0] == [] ) {
201: R = cdr(R); break;
202: }
203: }
1.5 noro 204: R = reverse(R);
1.3 noro 205: if ( DP ) return R;
206: else return map(dpmlisttollist,R,V);
207: }
208:
209: def dpmlisttollist(F,V)
210: {
211: return map(dpm_dtol,F,V);
212: }
213:
214: def sres(F,V,H,Ord)
215: {
216: if ( type(DP=getopt(dp)) == -1 ) DP = 0;
217: dpm_set_schreyer(0);
218: dp_ord(Ord);
219: K = length(F);
220: K = length(F);
221: for ( I = 0; I < K; I++ ) if ( F[I] ) break;
222: if ( I == K ) return [[],[],[]];
223: L = todpmlist(F,V);
224: F = L[0]; N = L[1];
1.6 ! noro 225: #if 0
1.3 noro 226: G = nd_gr(F,V,H,[0,Ord]|dp=1);
1.6 ! noro 227: #else
! 228: G = nd_gr_trace(F,V,H,1,[0,Ord]|dp=1);
! 229: #endif
1.3 noro 230: G = reverse(G);
231: R = [G];
232: dp_ord([0,Ord]);
233: while ( 1 ) {
234: S = dpm_schreyer_base(R[0]);
235: print(["length",length(S)]);
236: if ( S == [] ) break;
237: else R = cons(S,R);
238: }
239: dp_ord([0,0]);
1.5 noro 240: R = reverse(R);
1.3 noro 241: if ( DP ) return R;
242: else return map(dpmlisttollist,R,V);
243: }
244:
245: def minsres(F,V,H,Ord)
246: {
247: if ( type(DP=getopt(dp)) == -1 ) DP = 0;
248: R = sres(F,V,H,Ord|dp=1);
1.5 noro 249: R = ltov(reverse(R));
1.3 noro 250: M = length(R);
251: for ( I = 0; I < M; I++ ) R[I] = map(dpm_sort,R[I]);
252: R = vtol(R);
253: for ( T = R, U = []; length(T) >= 2; ) {
254: Z = dpm_simplify_syz(T[0],T[1]);
255: U = cons(Z[0],U);
256: T = cons(Z[1],cdr(cdr(T)));
257: }
258: U = cons(T[0],U);
259: if ( DP ) return U;
260: else return map(dpmlisttollist,U,V);
261: }
262:
263: def ordcheck(D,V)
264: {
265: B=map(dpm_dtol,D,V);
266: dp_ord([1,0]);
267: D = map(dpm_sort,D);
268: D1 = map(dpm_ltod,B,V);
269: if ( D != D1 )
270: print("afo");
271: return D1;
272: }
273:
274: def complsb(A,B)
275: {
276: AL = A[3]; AD = A[4];
277: BL = B[3]; BD = B[4];
278: if ( AD > BD ) return 1;
279: else if ( AD < BD ) return -1;
280: else if ( AL < BL ) return 1;
281: else if ( AL > BL ) return -1;
282: else return 0;
283: }
284:
285: def complsb_sd(A,B)
286: {
287: AL = A[3]; AD = A[4];
288: BL = B[3]; BD = B[4];
289: if ( AD-AL > BD-BL ) return 1;
290: else if ( AD-AL < BD-BL ) return -1;
291: else if ( AL > BL ) return 1;
292: else if ( AL < BL ) return -1;
293: else return 0;
294: }
295:
296: def sortlsb(A)
297: {
298: S = [];
299: N = length(A);
300: for ( I = 0; I < N; I++ ) S = append(cdr(vtol(A[I])),S);
301: // S = qsort(S,newsyz.complsb_sd);
302: S = qsort(S,newsyz.complsb);
303: return S;
304: }
305: /* D=[in,sum,deg,lev,ind] */
306: /* B=[0 B1 B2 ...] */
307: /* C=[0 C1 C2 ...] */
308: /* H=[0 H1 H2 ...] */
309: /* Z=[0 Z1 Z2 ...] */
310: /* Zi=[0 L1 L2 ...] */
311: /* Lj=[I1,I2,...] */
312: /* One=<<0,...,0>> */
313:
314: extern Rtime,Stime,Ptime$
315:
316: // B is sorted wrt positon
317: def find_pos(HF,B)
318: {
319: Len = length(B);
320: Pos = dpm_hp(HF);
321: First = 0; Last = Len-1;
322: if ( B[First][0] == HF ) return B[First][2];
323: if ( B[Last][0] == HF ) return B[Last][2];
324: while ( 1 ) {
325: Mid = idiv(First+Last,2);
326: PosM = dpm_hp(B[Mid][0]);
327: if ( PosM == Pos ) break;
328: else if ( PosM < Pos )
329: First = Mid;
330: else
331: Last = Mid;
332: }
333: for ( I = Mid; I < Len && dpm_hp(B[I][0]) == Pos; I++ )
334: if ( HF == B[I][0] ) return B[I][2];
335: for ( I = Mid-1; I >= 1 && dpm_hp(B[I][0]) == Pos; I-- )
336: if ( HF == B[I][0] ) return B[I][2];
337: error("find_pos : cannot happen");
338: }
339:
340: def reduce(D,B,Bpos,C,H,Z,K,Kind,G,One,Top)
341: {
342: M = D[0]; Ind = D[2]; I = D[3];
343: if ( I == 1 ) {
344: C[I][Ind] = G[Ind];
345: dpm_insert_to_zlist(Z[1],dpm_hp(G[Ind]),Ind);
346: H[I][Ind] = G[Ind];
347: } else {
348: /* M is in F(I-1) => phi(M) is in F(I-2) */
349: /* reduction in F(I-2) */
350: T0 = time()[0];
351: dpm_set_schreyer_level(I-2);
352: CI = C[I-1]; ZI = Z[I-1]; BI = B[I-1]; BposI = Bpos[I-1];
353: Len = size(CI);
354: T=dpm_hc(M); EM=dpm_hp(M);
355: XiM = T*dpm_ht(CI[EM]);
356: EXiM = dpm_hp(XiM);
357: ZIE = ZI[EXiM];
358: for ( ZIE = ZI[EXiM]; ZIE != []; ZIE = cdr(ZIE) ) {
359: J = car(ZIE);
360: if ( J > EM && dpm_redble(XiM,dpm_ht(CI[J])) ) break;
361: }
362: Ptime += time()[0]-T0;
363: T0 = time()[0];
364: QR = dpm_sp_nf(CI,ZI,EM,J|top=Top);
365: Rtime += time()[0]-T0;
366: G = QR[0]; F = QR[1];
367: if ( F ) {
368: HF = dpm_ht(F); EF = dpm_hp(HF);
369: /* find HF in B[I-1] */
370: T0 = time()[0];
371: J = find_pos(HF,BposI);
372: Stime += time()[0]-T0;
373: /* F=Ret[0]*Ret[1] */
374: Ret = dpm_remove_cont(F);
375: CI[J] = Ret[1];
376: Kind[I-1][J] = 2;
377: dpm_insert_to_zlist(ZI,EF,J);
378: dpm_set_schreyer_level(I-1);
379: Tail = -Ret[0]*dpm_dptodpm(One,J);
380: G += Tail;
381: dpm_set_schreyer_level(0);
382: Ret = dpm_remove_cont(G); G = Ret[1];
383: C[I][Ind] = G;
384: Kind[I][Ind] = 1;
385: K[I] = cons([G,Tail],K[I]);
386: dpm_insert_to_zlist(Z[I],EM,Ind);
387: /* level <- 0 */
388: BI[J][3] = 0;
389: } else {
390: Kind[I][Ind] = 0;
391: Ret = dpm_remove_cont(G); G = Ret[1];
392: H[I][Ind] = G;
393: C[I][Ind] = G;
394: dpm_insert_to_zlist(Z[I],EM,Ind);
395: }
396: }
397: }
398:
399: def lres_setup(F,V,H,Ord)
400: {
401: dpm_set_schreyer(0);
402: dp_ord(Ord);
403: K = length(F);
404: if ( type(F[0]) <= 2 ) {
405: // F is a polynimial list
406: F = map(dp_ptod,F,V);
407: F = map(dpm_dptodpm,F,1);
408: N = 1;
409: } else if ( type(F[0]) == 9 ) {
410: // F is a dpoly list
411: F = map(dpm_dptodpm,F,1);
412: N = 1;
413: } else if ( type(F[0]) == 4 ) {
414: // F is a list of poly lists
415: N = length(F[0]);
416: F = map(dpm_ltod,F,V);
417: } else if ( type(F[0]) == 26 ) {
418: // F is a DPM list
419: for ( N = 0, T = F; T != []; T = cdr(T) ) {
420: for ( A = car(T); A; A = dpm_rest(A) ) {
421: N1 = dpm_hp(A);
422: if ( N1 > N ) N = N1;
423: }
424: }
425: } else {
426: error("lres_setup: arugument type is invalid.");
427: }
1.6 ! noro 428: dp_ord([0,Ord]);
! 429: F = map(dpm_sort,F);
! 430: #if 0
1.3 noro 431: G = nd_gr(F,V,H,[0,Ord]|dp=1);
1.6 ! noro 432: #else
! 433: G = nd_gr_trace(F,V,H,1,[0,Ord]|dp=1);
! 434: #endif
1.3 noro 435: G = reverse(G);
436: dp_ord([0,Ord]);
437: One = dp_ptod(1,V);
438: return [G,One];
439: }
440:
441: def dpm_sort1(L)
442: {
443: return [dpm_sort(L[0]),L[1]];
444: }
445:
446: def comp_pos(A,B)
447: {
448: PA = dpm_hp(A[0]); PB = dpm_hp(B[0]);
449: if ( PA > PB ) return 1;
450: else if ( PA < PB ) return -1;
451: else return 0;
452: }
453:
454: def lres(F,V,H,Ord)
455: {
456: T0 = time();
457: if ( type(Top=getopt(top)) == -1 ) Top = 0;
1.5 noro 458: if ( type(DP=getopt(dp)) == -1 ) DP = 0;
1.3 noro 459: if ( type(NoSimpK=getopt(nosimpk)) == -1 ) NoSimpK = 0;
460: if ( type(NoPreProj=getopt(nopreproj)) == -1 ) NoPreProj = 0;
461: Rtime = Stime = Ptime = 0;
462: L = lres_setup(F,V,H,Ord);
463: G = L[0];
464: One = L[1];
465: F = dpm_schreyer_frame(G);
466: G = ltov(cons(0,L[0]));
467: F = reverse(F);
468: F = ltov(F);
469: N = length(F);
470: for ( I = 0; I < N; I++ ) {
471: FI = F[I]; FI[0] = [];
472: FI = map(ltov,FI);
473: F[I] = FI;
474: }
475: R = sortlsb(F);
476: B = vector(N+1);
477: Bpos = vector(N+1);
478: C = vector(N+1);
479: H = vector(N+1);
480: Z = vector(N+1);
481: K = vector(N+1);
482: L = vector(N+1);
483: K = vector(N+1);
484: D = vector(N+1);
485: Kind = vector(N+1);
486:
487: for ( I = 1; I <= N; I++ ) {
488: FI = F[I-1]; Len = length(FI);
489: B[I] = FI;
490: T = vector(Len-1);
491: for ( J = 1; J < Len; J++ ) T[J-1] = FI[J];
492: Bpos[I] = qsort(T,newsyz.comp_pos);
493: C[I] = vector(Len);
494: H[I] = vector(Len);
495: Kind[I] = vector(Len);
496: K[I] = [];
497: Max = 0;
498: for ( J = 1; J < Len; J++ )
499: if ( (Pos = dpm_hp(FI[J][0])) > Max ) Max = Pos;
500: Z[I] = ZI = vector(Max+1);
501: for ( J = 1; J <= Max; J++ ) ZI[J] = [];
502: }
503: T1 = time(); Ftime = T1[0]-T0[0];
504: R = ltov(R); Len = length(R);
505: print(["Len",Len]);
506: for ( I = 0, NF = 0; I < Len; I++ ) {
507: if ( !((I+1)%100) ) print(".",2);
508: if ( !((I+1)%10000) ) print(I+1);
509: if ( !R[I][3] ) continue;
510: NF++;
511: reduce(R[I],B,Bpos,C,H,Z,K,Kind,G,One,Top);
512: }
513: print("");
514: print(["NF",NF]);
515: T0 = time();
516: dpm_set_schreyer_level(0);
517: D[1] = map(dpm_sort,H[1]);
518: for ( I = 2; I <= N; I++ ) {
519: // HTT = [Head,Tab,TailTop]
520: HTT = create_base_ord(K[I],length(Kind[I-1]));
521: Head = HTT[0]; Tab = HTT[1]; TailTopPos = HTT[2];
522: if ( !NoPreProj )
523: Tab = map(remove_k,map(dpm_sort,Tab),Kind[I-1]);
524: else
525: Tab = map(dpm_sort,Tab);
526: TailTop = dpm_dptodpm(One,TailTopPos);
527: if ( !NoSimpK ) {
528: print("simplify_k "+rtostr(I)+"...",2);
529: simplify_k(Head,Tab,TailTop,One);
530: print("done");
531: }
532: HI = map(remove_k,map(dpm_sort,H[I]),Kind[I-1]);
533: Len = length(HI);
534: print("simplify_by_k "+rtostr(I)+"...",2);
535: D[I] = vector(Len);
536: for ( J = 0; J < Len; J++ ) {
537: D[I][J] = simplify_by_k(HI[J],Tab,TailTop,One);
538: if ( NoPreProj )
539: D[I][J] = remove_k(D[I][J],Kind[I-1]);
540: }
541: print("done");
542: }
543: dp_ord([1,0]);
544: T1 = time();
545: print(["Frame",Ftime,"Prep",Ptime,"Reduce",Rtime,"Search",Stime,"Minimalize",T1[0]-T0[0]]);
546: // return [C,H,K,Kind,D];
547: D = compress_h(D);
1.5 noro 548: if ( DP ) return D;
549: else return vtol(map(dpmlisttollist,D,V));
1.3 noro 550: }
551:
552: def create_base_ord(K,N)
553: {
554: Tab = vector(N);
555: Ks = [];
556: for ( T = K; T != []; T = cdr(T) ) {
557: Ks = cons(I=dpm_hp(car(T)[1]),Ks);
558: Tab[I] = car(T)[0];
559: }
560: Others = [];
561: for ( I = N-1; I >= 1; I-- )
562: if ( !Tab[I] ) Others = cons(I,Others);
563: Ks = reverse(Ks);
564: dp_ord([4,append(Ks,Others),0]);
565: return [ltov(Ks),Tab,car(Others)];
566: }
567:
568: /* Head = [i1,i2,...], Tab[i1] = K[0], Tab[i2] = K[1],... */
569:
570: def simplify_k(Head,Tab,TailTop,One)
571: {
572: N = length(Tab);
573: Len = length(Head);
574: for ( I = Len-1; I >= 0; I-- ) {
575: M = 1;
576: R = Tab[Head[I]];
577: H = dpm_hm(R); R = dpm_rest(R);
578: while ( R && dpm_dptodpm(One,dpm_hp(R)) > TailTop ) {
579: Pos = dpm_hp(R); Red = Tab[Pos];
580: CRed = dp_hc(dpm_hc(Red));
581: CR = dpm_extract(R,Pos);
582: L = dpm_remove_cont(CRed*R-CR*Red);
583: R = L[1]; M *= CRed/L[0];
584: }
585: Tab[Head[I]] = M*H+R;
586: }
587: }
588:
589: def simplify_by_k(F,Tab,TailTop,One)
590: {
591: R = F;
592: M = 1;
593: while ( R && dpm_dptodpm(One,dpm_hp(R)) > TailTop ) {
594: Pos = dpm_hp(R); Red = Tab[Pos];
595: CRed = dp_hc(dpm_hc(Red));
596: CR = dpm_extract(R,Pos);
597: L = dpm_remove_cont(CRed*R-CR*Red);
598: R = L[1]; M *= CRed/L[0];
599: }
600: return (1/M)*R;
601: }
602:
603: /* Kind[I]=0 => phi(ei)\in H, Kind[I]=1 => phi(ei)\in K */
604: def remove_k(F,Kind)
605: {
606: R = [];
607: for ( T = F; T; T = dpm_rest(T) )
608: if ( Kind[dpm_hp(T)] != 1 ) R = cons(dpm_hm(T),R);
609: for ( S = 0; R != []; R = cdr(R) )
610: S += car(R);
611: return S;
612: }
613:
614: def remove_k1(F,Kind)
615: {
616: return [remove_k(F[0],Kind),F[1]];
617: }
618:
619: def extract_nonzero(A)
620: {
621: N = length(A);
622: R = [];
623: for ( C = I = 0; I < N; I++ )
624: if ( A[I] ) R=cons(A[I],R);
625: return reverse(R);;
626: }
627:
628: def nonzero(A)
629: {
630: N = length(A);
631: for ( C = I = 0; I < N; I++ )
632: if ( A[I] ) C++;
633: return C;
634: }
635:
636: def phi(C,F)
637: {
638: R = 0;
639: for ( T = F; T; T = dpm_rest(T) ) {
640: Coef = dpm_hc(T); Pos = dpm_hp(T);
1.6 ! noro 641: R += Coef*C[Pos-1];
1.3 noro 642: }
643: return R;
1.1 noro 644: }
645:
1.3 noro 646: def syz_check(H,I)
1.1 noro 647: {
1.3 noro 648: HI = H[I];
649: // dpm_set_schreyer_level(I-1);
650: HI1 = H[I+1];
651: Len = size(HI1)[0];
652: for ( J = 1; J < Len; J++ ) {
653: F = HI1[J];
654: R = phi(HI,F);
655: if ( R ) print([J,"NG"]);
656: }
1.1 noro 657: }
658:
1.3 noro 659: // for compressed C
660: def phi0(C,F)
1.1 noro 661: {
1.3 noro 662: R = 0;
663: for ( T = F; T; T = dpm_rest(T) ) {
664: Coef = dpm_hc(T); Pos = dpm_hp(T);
665: R += Coef*C[Pos-1];
666: }
667: return R;
1.1 noro 668: }
669:
1.3 noro 670: // for compressed H
671: def syz_check0(H,I)
1.1 noro 672: {
1.3 noro 673: HI = H[I];
674: // dpm_set_schreyer_level(I-1);
675: HI1 = H[I+1];
676: Len = length(HI1);
677: for ( J = 0; J < Len; J++ ) {
678: F = HI1[J];
679: R = phi0(HI,F);
680: if ( R ) print([J,"NG"]);
681: }
1.1 noro 682: }
683:
1.3 noro 684: // renumber position
685: def renumber_pos(F,Tab)
1.1 noro 686: {
1.3 noro 687: L = [];
688: for ( T = F; T; T = dpm_rest(T) )
689: L = cons(dpm_hm(T),L);
690: R = 0;
691: for ( T = L; T != []; T = cdr(T) )
692: R += dpm_dptodpm(dpm_hc(car(T)),Tab[dpm_hp(car(T))]);
693: return R;
1.1 noro 694: }
695:
1.3 noro 696: // compress H1 and renumber H
697: def compress(H,H1)
1.1 noro 698: {
1.3 noro 699: // create index table for H1
700: L1 = length(H1);
701: Tmp = vector(L1);
702: Tab = vector(L1);
703: for ( I = 1, J = 1; I < L1; I++ )
704: if ( H1[I] ) {
705: Tab[I] = J;
706: Tmp[J++] = H1[I];
707: }
708: NH1 = vector(J);
709: for ( I = 1; I < J; I++ )
710: NH1[I] = Tmp[I];
711: if ( H )
712: H = map(renumber_pos,H,Tab);
713: return [H,NH1];
1.1 noro 714: }
715:
1.3 noro 716: def compress_h(H)
1.1 noro 717: {
1.3 noro 718: // H = [0,H[1],...,H[L-1]]
719: L = length(H);
720: NH = vector(L-1);
721: H1 = H[1];
722: for ( I = 2; I < L; I++ ) {
723: R = compress(H[I],H1);
724: H1 = R[0];
725: NH[I-2] = cdr(vtol(R[1]));
726: }
727: R = compress(0,H1);
728: NH[L-2] = cdr(vtol(R[1]));
729: return NH;
1.1 noro 730: }
1.3 noro 731: endmodule$
1.1 noro 732: end$
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>