![[BACK]](/icons/cvsweb/back.gif){kind=icon}
Return to de.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 de.rr CVS log | Up to [local] / OpenXM / src / asir-contrib / testing / noro |
1.1 noro 1: module de;
2:
3: localf split_lexgb;
4: localf sort_lex_dec,sort_lex_inc;
5: localf inverse_or_split, linear_dim;
6: localf dp_monic_mod,monic_gb;
7: localf membership_test;
8: localf dp_chrem,intdptoratdp,intdpltoratdpl;
9: localf comp_by_ht,dp_gr_mod,gr_chrem;
10: localf construct_sqfrbasis;
1.2 ! noro 11: localf split;
1.1 noro 12:
13: /*
14: * G : a 0-dim lex gb, reduced
15: * if there exists a non-monic g = g(x')x^k+...
16: * then g(x') is a zero divisor and Id(G) splits
17: */
18:
19: def split_lexgb(G,V)
20: {
21: G = map(ptozp,G);
22: G = sort_lex_dec(G,V);
23: TG = G;
24: for ( ; TG != []; TG = cdr(TG) ) {
25: F = car(TG);
26: for ( TV = V; !deg(F,car(TV)); TV = cdr(TV) );
27: VF = car(TV);
28: DF = deg(F,VF);
29: CF = coef(F,DF,VF);
30: if ( type(CF) == 2 ) {
31: L = inverse_or_split(V,G,CF);
32: R = map(split_lexgb,L,V);
33: return append(R[0],R[1]);
34: }
35: }
36: return [G];
37: }
38:
39: /*
40: * V : a variable list
41: * Id : a 0-dim radical ideal represented by a lex basis
42: * F : a poly
43: * if F is a unit of Q[V]/Id, then 1/F is returned
44: * else F is a zero divisor and Id = (Id:F)\cap(Id+<F>)
45: * [GB(Id:F),GB(Id+<F>)] is returned.
46: */
47:
48: def inverse_or_split(V,Id,F)
49: {
50: Id = map(ptozp,Id);
51: N = length(V);
52: dp_ord(2);
53: set_field(Id,V,2);
54: DF = dptodalg(dp_ptod(F,V));
55: Ret = inv_or_split_dalg(DF);
56: if ( type(Ret) == 1 ) {
57: /* Ret = 1/F */
58: Dp = dalgtodp(Ret);
59: return dp_dtop(Dp[0],V)/Dp[1];
60: } else {
61: /* Ret = GB(Id:F) */
62: /* compute GB(Id+<f>) */
63: Gquo = append(map(ptozp,map(dp_dtop,Ret,V)),Id);
64: /* inter-reduction */
65: Gquo = nd_gr_postproc(Gquo,V,0,2,0);
66: DTotal = linear_dim(Id,V,2);
67: Dquo = linear_dim(Gquo,V,2);
68: Drem = DTotal-Dquo;
69: B = cons(F,Id);
70: Grem = gr_chrem(B,V,2|dim=Drem);
71: return [map(ptozp,Gquo),map(ptozp,Grem)];
72: }
73: }
74:
1.2 ! noro 75: def split(V,Id,F,Ord)
! 76: {
! 77: Id = map(ptozp,Id);
! 78: N = length(V);
! 79: dp_ord(Ord);
! 80: set_field(Id,V,Ord);
! 81: DF = dptodalg(dp_ptod(F,V));
! 82: Ret = inv_or_split_dalg(DF);
! 83: /* Ret = GB(Id:F) */
! 84: /* compute GB(Id+<f>) */
! 85: Gquo = append(map(ptozp,map(dp_dtop,Ret,V)),Id);
! 86: /* inter-reduction */
! 87: Gquo = nd_gr_postproc(Gquo,V,0,Ord,0);
! 88: DTotal = linear_dim(Id,V,Ord);
! 89: Dquo = linear_dim(Gquo,V,Ord);
! 90: Drem = DTotal-Dquo;
! 91: B = cons(F,Id);
! 92: Grem = gr_chrem(B,V,Ord|dim=Drem);
! 93: return [map(ptozp,Gquo),map(ptozp,Grem)];
! 94: }
! 95:
1.1 noro 96: /* add F(X,V) to Id(B) */
97: /* returns a list of splitted ideals */
98: /* B should be a triangular basis */
99:
100: def construct_sqfrbasis(F,X,B,V)
101: {
102: if ( type(F) == 1 )
103: return [];
104: B = sort_lex_dec(B,V);
105: V1 = cons(X,V);
106: F = nd_nf(F,reverse(B),cons(X,V),2,0);
107: D = deg(F,X);
108: H = coef(F,D,X);
109: if ( type(H) == 2 ) {
110: Ret = inverse_or_split(V,B,H);
111: if ( type(Ret) == 4 ) {
112: /* H != 0 on Id_nz, H = 0 on Id_z */
113: B0=construct_sqfrbasis(F,X,Ret[0],V);
114: B1=construct_sqfrbasis(F,X,Ret[1],V);
115: return append(B0,B1);
116: } else
117: F = nd_nf(F*Ret,reverse(B),cons(X,V),2,0);
118: }
119: B1 = cons(F,B);
120: /* F is monic */
121: M = minipoly(B1,V1,2,X,zzz);
122: S = sqfr(M); S = cdr(S);
123: if ( length(S) == 1 && car(S)[1] == 1 )
124: return [cons(F,B)];
125: else {
126: R = [];
127: for ( T = S; T != []; T = cdr(T) ) {
128: G = nd_gr_trace(cons(subst(car(T)[0],zzz,X),B1),V1,1,1,2);
129: R1 = split_lexgb(G,V1);
130: R = append(R1,R);
131: }
132: return R;
133: }
134: }
135:
136: def sort_lex_dec(B,V)
137: {
138: dp_ord(2);
139: B = map(dp_ptod,B,V);
140: B = vtol(qsort(ltov(B),comp_by_ht));
141: B = map(dp_dtop,B,V);
142: return reverse(B);
143: }
144:
145: def sort_lex_inc(B,V)
146: {
147: dp_ord(2);
148: B = map(dp_ptod,B,V);
149: B = vtol(qsort(ltov(B),comp_by_ht));
150: B = map(dp_dtop,B,V);
151: return B;
152: }
153:
154: def linear_dim(G,V,Ord)
155: {
156: dp_ord(Ord);
157: MB = dp_mbase(map(dp_ptod,G,V));
158: return length(MB);
159: }
160:
161: def membership_test(B,G,V,O)
162: {
163: B = map(ptozp,B);
164: G = map(ptozp,G);
165: for ( T = B; T != []; T = cdr(T) )
166: if ( nd_nf(car(T),G,V,O,0) ) return 0;
167: return 1;
168: }
169:
170: def gr_chrem(B,V,O)
171: {
172: /* linear dimension */
173: Dim = getopt(dim);
174: if ( type(Dim) == -1 ) Dim = -1;
175:
176: B = map(ptozp,B);
177: G = []; HS = []; Mod = 1;
178: for ( I = 0; ; I++ ) {
179: P = lprime(I);
180: GM = nd_gr(B,V,P,O);
181: if ( Dim >= 0 )
182: if ( linear_dim(GM,V,O) != Dim ) continue;
183: L = monic_gb(GM,V,O,P); GM = L[0]; HSM = L[1];
184: G = dp_chrem(G,HS,Mod,GM,HSM,P);
185: Mod *= P;
186: if ( G != [] )
187: HS = HSM;
188: M = idiv(isqrt(2*Mod),2);
189: R1 = intdpltoratdpl(G,Mod,M);
190: if ( R1 ) {
191: if ( Found && R == R1
192: && (GB=nd_gr_postproc(map(dp_dtop,R,V),V,0,O,1))
193: && membership_test(B,GB,V,O) )
194: break;
195: else {
196: R = R1; Found = 1;
197: }
198: }
199: }
200: return GB;
201: }
202:
203: def comp_by_ht(A,B)
204: {
205: HA = dp_ht(A); HB = dp_ht(B);
206: if ( HA > HB )
207: return 1;
208: else if ( HA < HB )
209: return -1;
210: else
211: return 0;
212: }
213:
214: def intdpltoratdpl(G,Mod,M)
215: {
216: for ( R = []; G != []; G = cdr(G) ) {
217: T = intdptoratdp(car(G),Mod,M);
218: if ( !T )
219: return 0;
220: else
221: R = cons(T,R);
222: }
223: R = reverse(R);
224: return vtol(qsort(newvect(length(R),R),comp_by_ht));
225: }
226:
227: def intdptoratdp(F,Mod,M)
228: {
229: for ( T = F, N = 0; T; T = dp_rest(T), N++ );
230: C = newvect(N);
231: for ( I = 0, T = F; I < N; T = dp_rest(T), I++ ) {
232: L = inttorat(dp_hc(T),Mod,M);
233: if ( !L )
234: return 0;
235: else
236: C[I] = (L[0]/L[1])*dp_ht(T);
237: }
238: for ( R = 0, I = N-1; I >= 0; I-- )
239: R += C[I];
240: return R;
241: }
242:
243: def dp_chrem(G,HS,Mod,GM,HSM,P)
244: {
245: if ( G == [] )
246: return GM;
247: if ( HS != HSM )
248: return [];
249: R = [];
250: M1 = inv(Mod,P);
251: ModM1 = Mod*M1;
252: for ( ; G != []; G = cdr(G), GM = cdr(GM) ) {
253: E = car(G); EM = car(GM);
254: E1 = E+(EM-E)*ModM1;
255: R = cons(E1,R);
256: }
257: return reverse(R);
258: }
259:
260: def monic_gb(G,V,O,P)
261: {
262: dp_ord(O); setmod(P);
263: D = map(dp_ptod,G,V);
264: D = map(dp_monic_mod,D,P);
265: D = vtol(qsort(newvect(length(D),D),comp_by_ht));
266: return [D,map(dp_ht,D)];
267: }
268:
269: def dp_monic_mod(F,P)
270: {
271: FP = dp_mod(F,P,[]);
272: return dp_rat(FP/dp_hc(FP));
273: }
274:
275: endmodule;
276: end$