1 | /**************************************** |
---|
2 | * Computer Algebra System SINGULAR * |
---|
3 | ****************************************/ |
---|
4 | /* $Id: syz2.cc,v 1.7 1999-10-22 11:14:19 obachman Exp $ */ |
---|
5 | /* |
---|
6 | * ABSTRACT: resolutions |
---|
7 | */ |
---|
8 | #include <limits.h> |
---|
9 | |
---|
10 | #include "mod2.h" |
---|
11 | #include "tok.h" |
---|
12 | #include "mmemory.h" |
---|
13 | #include "syz.h" |
---|
14 | #include "polys.h" |
---|
15 | #include "febase.h" |
---|
16 | #include "kstd1.h" |
---|
17 | #include "kutil.h" |
---|
18 | #include "stairc.h" |
---|
19 | #include "ipid.h" |
---|
20 | #include "cntrlc.h" |
---|
21 | #include "ipid.h" |
---|
22 | #include "intvec.h" |
---|
23 | #include "ipshell.h" |
---|
24 | #include "limits.h" |
---|
25 | #include "numbers.h" |
---|
26 | #include "modulop.h" |
---|
27 | #include "ideals.h" |
---|
28 | #include "intvec.h" |
---|
29 | #include "ring.h" |
---|
30 | #include "lists.h" |
---|
31 | #include "kbuckets.h" |
---|
32 | #include "polys-comp.h" |
---|
33 | |
---|
34 | //#define SHOW_PROT |
---|
35 | //#define SHOW_RED |
---|
36 | //#define SHOW_HILB |
---|
37 | //#define SHOW_RESULT |
---|
38 | //#define INVERT_PAIRS |
---|
39 | //#define SHOW_CRIT |
---|
40 | //#define SHOW_SPRFL |
---|
41 | #define USE_CHAINCRIT |
---|
42 | /*--- some heuristics to optimize the pair sets---*/ |
---|
43 | /*--- chose only one (!!!) at the same time ------*/ |
---|
44 | //#define USE_HEURISTIC1 |
---|
45 | #define USE_HEURISTIC2 |
---|
46 | |
---|
47 | extern void rSetmS(poly p, int* Components, long* ShiftedComponents); |
---|
48 | #ifdef SHOW_CRIT |
---|
49 | static int crit; |
---|
50 | static int crit1; |
---|
51 | static int spfl; |
---|
52 | static int cons_pairs; |
---|
53 | static int crit_fails; |
---|
54 | #endif |
---|
55 | typedef struct sopen_pairs open_pairs; |
---|
56 | typedef open_pairs* crit_pairs; |
---|
57 | struct sopen_pairs |
---|
58 | { |
---|
59 | crit_pairs next; |
---|
60 | int first_poly; |
---|
61 | int second_poly; |
---|
62 | }; |
---|
63 | /*3 |
---|
64 | * computes pairs from the new elements (beginning with the element newEl) |
---|
65 | * in the module index |
---|
66 | */ |
---|
67 | static void syCreateNewPairs_Hilb(syStrategy syzstr, int index, |
---|
68 | int actdeg) |
---|
69 | { |
---|
70 | SSet temp; |
---|
71 | SObject tso; |
---|
72 | poly toHandle,tsyz=NULL,p,pp; |
---|
73 | int r1,r2=0,rr,l=(*syzstr->Tl)[index]; |
---|
74 | int i,j,r=0,ti; |
---|
75 | BOOLEAN toComp=FALSE; |
---|
76 | #ifdef USE_CHAINCRIT |
---|
77 | crit_pairs cp=NULL,tcp; |
---|
78 | #endif |
---|
79 | actdeg += index; |
---|
80 | long * ShiftedComponents = syzstr->ShiftedComponents[index-1]; |
---|
81 | int* Components = syzstr->truecomponents[index-1]; |
---|
82 | |
---|
83 | while ((l>0) && ((syzstr->resPairs[index])[l-1].lcm==NULL)) l--; |
---|
84 | rr = l-1; |
---|
85 | while ((rr>=0) && (((syzstr->resPairs[index])[rr].p==NULL) || |
---|
86 | ((syzstr->resPairs[index])[rr].order>actdeg))) rr--; |
---|
87 | r2 = rr+1; |
---|
88 | while ((rr>=0) && ((syzstr->resPairs[index])[rr].order==actdeg) |
---|
89 | && ((syzstr->resPairs[index])[rr].syzind<0)) |
---|
90 | { |
---|
91 | rr--; |
---|
92 | r++; |
---|
93 | } |
---|
94 | if (r==0) return; |
---|
95 | ideal nP=idInit(l,syzstr->res[index]->rank); |
---|
96 | #ifdef INVERT_PAIRS |
---|
97 | r1 = rr+1; |
---|
98 | #else |
---|
99 | r1 = r2-1; |
---|
100 | #endif |
---|
101 | /*---------- there are new pairs ------------------------------*/ |
---|
102 | loop |
---|
103 | { |
---|
104 | /*--- chose first new elements --------------------------------*/ |
---|
105 | toComp = FALSE; |
---|
106 | toHandle = (syzstr->resPairs[index])[r1].p; |
---|
107 | if (toHandle!=NULL) |
---|
108 | { |
---|
109 | int tc=pGetComp(toHandle); |
---|
110 | (syzstr->resPairs[index])[r1].syzind = 0; |
---|
111 | for (i=0; i<r1;i++) |
---|
112 | { |
---|
113 | if (((syzstr->resPairs[index])[i].p!=NULL) && |
---|
114 | (pGetComp((syzstr->resPairs[index])[i].p)==tc)) |
---|
115 | { |
---|
116 | #ifdef USE_CHAINCRIT |
---|
117 | tcp = cp; |
---|
118 | if (tcp!=NULL) |
---|
119 | { |
---|
120 | while ((tcp!=NULL) && |
---|
121 | ((tcp->first_poly!=i)||(tcp->second_poly!=r1))) tcp = tcp->next; |
---|
122 | } |
---|
123 | if (tcp==NULL) |
---|
124 | { |
---|
125 | #endif |
---|
126 | p = pOne(); |
---|
127 | pLcm((syzstr->resPairs[index])[i].p,toHandle,p); |
---|
128 | rSetmS(p,Components, ShiftedComponents); |
---|
129 | j = 0; |
---|
130 | while (j<i) |
---|
131 | { |
---|
132 | if (nP->m[j]!=NULL) |
---|
133 | { |
---|
134 | if (pDivisibleBy2(nP->m[j],p)) |
---|
135 | { |
---|
136 | pDelete(&p); |
---|
137 | p = NULL; |
---|
138 | break; |
---|
139 | } |
---|
140 | else if (pDivisibleBy2(p,nP->m[j])) |
---|
141 | { |
---|
142 | pDelete(&(nP->m[j])); |
---|
143 | nP->m[j] = NULL; |
---|
144 | } |
---|
145 | #ifdef USE_CHAINCRIT |
---|
146 | else |
---|
147 | { |
---|
148 | poly p1,p2; |
---|
149 | int ip=pVariables; |
---|
150 | p1 = pDivide(p,(syzstr->resPairs[index])[i].p); |
---|
151 | p2 = pDivide(nP->m[j],(syzstr->resPairs[index])[j].p); |
---|
152 | while ((ip>0) && (pGetExp(p1,ip)*pGetExp(p2,ip)==0)) ip--; |
---|
153 | if (ip==0) |
---|
154 | { |
---|
155 | #ifdef SHOW_SPRFL |
---|
156 | Print("Hier: %d, %d\n",j,i); |
---|
157 | #endif |
---|
158 | if (i>rr) |
---|
159 | { |
---|
160 | tcp=(crit_pairs)Alloc0(sizeof(sopen_pairs)); |
---|
161 | tcp->next = cp; |
---|
162 | tcp->first_poly = j; |
---|
163 | tcp->second_poly = i; |
---|
164 | cp = tcp; |
---|
165 | tcp = NULL; |
---|
166 | } |
---|
167 | else |
---|
168 | { |
---|
169 | ti=0; |
---|
170 | while ((ti<l) && (((syzstr->resPairs[index])[ti].ind1!=j)|| |
---|
171 | ((syzstr->resPairs[index])[ti].ind2!=i))) ti++; |
---|
172 | if (ti<l) |
---|
173 | { |
---|
174 | #ifdef SHOW_SPRFL |
---|
175 | Print("gefunden in Mod %d: ",index); pWrite((syzstr->resPairs[index])[ti].lcm); |
---|
176 | #endif |
---|
177 | syDeletePair(&(syzstr->resPairs[index])[ti]); |
---|
178 | #ifdef SHOW_CRIT |
---|
179 | crit1++; |
---|
180 | #endif |
---|
181 | toComp = TRUE; |
---|
182 | } |
---|
183 | } |
---|
184 | } |
---|
185 | pDelete(&p1); |
---|
186 | pDelete(&p2); |
---|
187 | } |
---|
188 | #endif |
---|
189 | } |
---|
190 | j++; |
---|
191 | } |
---|
192 | if (p!=NULL) |
---|
193 | { |
---|
194 | nP->m[i] = p; |
---|
195 | } |
---|
196 | #ifdef USE_CHAINCRIT |
---|
197 | } |
---|
198 | else |
---|
199 | { |
---|
200 | #ifdef SHOW_CRIT |
---|
201 | crit1++; |
---|
202 | #endif |
---|
203 | } |
---|
204 | #endif |
---|
205 | } |
---|
206 | } |
---|
207 | if (toComp) syCompactify1(syzstr->resPairs[index],&l,r1); |
---|
208 | for (i=0;i<r1;i++) |
---|
209 | { |
---|
210 | if (nP->m[i]!=NULL) |
---|
211 | { |
---|
212 | tso.lcm = p = nP->m[i]; |
---|
213 | nP->m[i] = NULL; |
---|
214 | tso.order = pTotaldegree(p); |
---|
215 | if ((syzstr->cw!=NULL) && (index>0) && (pGetComp(p)>0)) |
---|
216 | { |
---|
217 | int ii=index-1,jj=pGetComp(p); |
---|
218 | while (ii>0) |
---|
219 | { |
---|
220 | jj = pGetComp(syzstr->res[ii]->m[jj-1]); |
---|
221 | ii--; |
---|
222 | } |
---|
223 | tso.order += (*syzstr->cw)[jj-1]; |
---|
224 | } |
---|
225 | tso.p1 = (syzstr->resPairs[index])[i].p; |
---|
226 | tso.p2 = toHandle; |
---|
227 | tso.ind1 = i; |
---|
228 | tso.ind2 = r1; |
---|
229 | tso.syzind = -1; |
---|
230 | tso.isNotMinimal = (poly)1; |
---|
231 | tso.p = NULL; |
---|
232 | tso.length = -1; |
---|
233 | currcomponents = syzstr->truecomponents[index]; |
---|
234 | currShiftedComponents = syzstr->ShiftedComponents[index]; |
---|
235 | rChangeSComps(currcomponents, |
---|
236 | currShiftedComponents, |
---|
237 | IDELEMS(syzstr->res[index])); // actueller index |
---|
238 | number coefgcd = |
---|
239 | nGcd(pGetCoeff(tso.p1),pGetCoeff(tso.p2)); |
---|
240 | tso.syz = pCopy((syzstr->resPairs[index])[i].syz); |
---|
241 | poly tt = pDivide(tso.lcm,tso.p1); |
---|
242 | pSetCoeff(tt,nDiv(pGetCoeff(tso.p1),coefgcd)); |
---|
243 | tso.syz = pMultT(tso.syz,tt); |
---|
244 | pDelete(&tt); |
---|
245 | coefgcd = nNeg(coefgcd); |
---|
246 | pp = pCopy((syzstr->resPairs[index])[r1].syz); |
---|
247 | tt = pDivide(tso.lcm,tso.p2); |
---|
248 | pSetCoeff(tt,nDiv(pGetCoeff(tso.p2),coefgcd)); |
---|
249 | pp = pMultT(pp,tt); |
---|
250 | pDelete(&tt); |
---|
251 | tso.syz = pAdd(pp,tso.syz); |
---|
252 | nDelete(&coefgcd); |
---|
253 | currcomponents = syzstr->truecomponents[index-1]; |
---|
254 | currShiftedComponents = syzstr->ShiftedComponents[index-1]; |
---|
255 | pSetComp(tso.lcm,pGetComp((syzstr->resPairs[index])[r1].syz)); |
---|
256 | rChangeSComps(currcomponents, |
---|
257 | currShiftedComponents, |
---|
258 | IDELEMS(syzstr->res[index])); // actueller index |
---|
259 | #ifdef SHOW_PROT |
---|
260 | Print("erzeuge Paar im Modul %d,%d mit: \n",index,tso.order); |
---|
261 | Print("poly1: ");pWrite(tso.p1); |
---|
262 | Print("poly2: ");pWrite(tso.p2); |
---|
263 | Print("syz: ");pWrite(tso.syz); |
---|
264 | Print("sPoly: ");pWrite(tso.p); |
---|
265 | PrintLn(); |
---|
266 | #endif |
---|
267 | syEnterPair(syzstr,&tso,&l,index); |
---|
268 | } |
---|
269 | } |
---|
270 | } |
---|
271 | #ifdef INVERT_PAIRS |
---|
272 | r1++; |
---|
273 | if (r1>=r2) break; |
---|
274 | #else |
---|
275 | r1--; |
---|
276 | if (r1<=rr) break; |
---|
277 | #endif |
---|
278 | } |
---|
279 | #ifdef USE_CHAINCRIT |
---|
280 | while (cp!=NULL) |
---|
281 | { |
---|
282 | tcp = cp; |
---|
283 | cp = cp->next; |
---|
284 | Free((ADDRESS)tcp,sizeof(sopen_pairs)); |
---|
285 | } |
---|
286 | #endif |
---|
287 | } |
---|
288 | |
---|
289 | /*3 |
---|
290 | * determines the place of a polynomial in the right ordered resolution |
---|
291 | * set the vectors of truecomponents |
---|
292 | */ |
---|
293 | static BOOLEAN syOrder_Hilb(poly p,syStrategy syzstr,int index, |
---|
294 | int realcomp) |
---|
295 | { |
---|
296 | int i=IDELEMS(syzstr->res[index-1])+1,j=0,k,tc,orc,ie; |
---|
297 | int *trind1=syzstr->truecomponents[index-1]; |
---|
298 | int *trind=syzstr->truecomponents[index]; |
---|
299 | long *shind=syzstr->ShiftedComponents[index]; |
---|
300 | poly pp; |
---|
301 | polyset o_r=syzstr->orderedRes[index]->m; |
---|
302 | BOOLEAN ret = FALSE; |
---|
303 | |
---|
304 | // if != 0, then new element can go into same component |
---|
305 | // i.e., we do not need to leave space in shifted components |
---|
306 | long same_comp = 0; |
---|
307 | |
---|
308 | if (p==NULL) return FALSE; |
---|
309 | if (realcomp==0) realcomp=1; |
---|
310 | ie = IDELEMS(syzstr->orderedRes[index]); |
---|
311 | while ((ie>0) && (syzstr->orderedRes[index]->m[ie-1]==NULL)) ie--; |
---|
312 | if (ie==0) |
---|
313 | { |
---|
314 | j = 0; |
---|
315 | } |
---|
316 | else |
---|
317 | { |
---|
318 | if (index>1) |
---|
319 | tc = trind1[pGetComp(p)]; |
---|
320 | else |
---|
321 | tc = pGetComp(p); |
---|
322 | loop |
---|
323 | { |
---|
324 | orc = pGetComp(o_r[j]); |
---|
325 | if (trind1[orc]>tc+1) break; |
---|
326 | else if (trind1[orc] == tc+1) |
---|
327 | { |
---|
328 | same_comp = 1; |
---|
329 | } |
---|
330 | j++; |
---|
331 | if (j==ie) break; |
---|
332 | } |
---|
333 | } |
---|
334 | if (j == ie) |
---|
335 | { |
---|
336 | // new element is the last in ordered module |
---|
337 | if (same_comp == 0) |
---|
338 | same_comp = SYZ_SHIFT_BASE; |
---|
339 | |
---|
340 | // test wheter we have enough space for new shifted component |
---|
341 | if ((LONG_MAX - same_comp) <= shind[ie]) |
---|
342 | { |
---|
343 | long new_space = syReorderShiftedComponents(shind, ie+1); |
---|
344 | assume((LONG_MAX - same_comp) > shind[ie]); |
---|
345 | ret = TRUE; |
---|
346 | if (TEST_OPT_PROT) Print("(T%u)", new_space); |
---|
347 | } |
---|
348 | |
---|
349 | // yes, then set new shifted component |
---|
350 | assume(ie == 0 || shind[ie] > 0); |
---|
351 | shind[ie+1] = shind[ie] + same_comp; |
---|
352 | } |
---|
353 | else |
---|
354 | { |
---|
355 | // new element must come in between |
---|
356 | // i.e. at place j+1 |
---|
357 | long prev, next; |
---|
358 | |
---|
359 | // test whether new component can get shifted value |
---|
360 | prev = shind[j]; |
---|
361 | next = shind[j+1]; |
---|
362 | assume(next > prev); |
---|
363 | if ((same_comp && prev + 2 >= next) || (!same_comp && next - prev < 4)) |
---|
364 | { |
---|
365 | long new_space = syReorderShiftedComponents(shind, ie+1); |
---|
366 | prev = shind[j]; |
---|
367 | next = shind[j+1]; |
---|
368 | assume((same_comp && prev + 2 < next) || (!same_comp && next - prev >= 4)); |
---|
369 | ret = TRUE; |
---|
370 | if (TEST_OPT_PROT) Print("(B%u)", new_space); |
---|
371 | } |
---|
372 | |
---|
373 | // make room for insertion of j+1 shifted component |
---|
374 | for (k=ie+1; k > j+1; k--) shind[k] = shind[k-1]; |
---|
375 | |
---|
376 | if (same_comp) |
---|
377 | { |
---|
378 | // can simply add one |
---|
379 | shind[j+1] = prev + 1; |
---|
380 | assume(shind[j+1] + 1 < shind[j+2]); |
---|
381 | } |
---|
382 | else |
---|
383 | { |
---|
384 | // need to leave more breathing room - i.e. value goes in |
---|
385 | // between |
---|
386 | shind[j+1] = prev + ((next - prev) >> 1); |
---|
387 | assume (shind[j] + 1 < shind[j+1] && shind[j+1] + 1 < shind[j+2]); |
---|
388 | } |
---|
389 | } |
---|
390 | |
---|
391 | if (o_r[j]!=NULL) |
---|
392 | { |
---|
393 | for (k=ie;k>j;k--) |
---|
394 | { |
---|
395 | o_r[k] = o_r[k-1]; |
---|
396 | } |
---|
397 | } |
---|
398 | o_r[j] = p; |
---|
399 | for (k=0;k<IDELEMS((syzstr->res)[index]);k++) |
---|
400 | { |
---|
401 | if (trind[k]>j) |
---|
402 | trind[k] += 1; |
---|
403 | } |
---|
404 | for (k=IDELEMS((syzstr->res)[index])-1;k>realcomp;k--) |
---|
405 | trind[k] = trind[k-1]; |
---|
406 | trind[realcomp] = j+1; |
---|
407 | return ret; |
---|
408 | } |
---|
409 | |
---|
410 | static void syHalfPair(poly syz, int newEl, syStrategy syzstr, int index) |
---|
411 | { |
---|
412 | SObject tso; |
---|
413 | int l=(*syzstr->Tl)[index]; |
---|
414 | |
---|
415 | while ((l>0) && ((syzstr->resPairs[index])[l-1].syz==NULL)) l--; |
---|
416 | if ((syzstr->cw!=NULL) && (index>0) && (pGetComp(syz)>0)) |
---|
417 | { |
---|
418 | int ii=index-1,jj=pGetComp(syz); |
---|
419 | while (ii>0) |
---|
420 | { |
---|
421 | jj = pGetComp(syzstr->res[ii]->m[jj-1]); |
---|
422 | ii--; |
---|
423 | } |
---|
424 | tso.order += (*syzstr->cw)[jj-1]; |
---|
425 | } |
---|
426 | tso.p1 = NULL; |
---|
427 | tso.p2 = NULL; |
---|
428 | tso.ind1 = 0; |
---|
429 | tso.ind2 = 0; |
---|
430 | tso.syzind = -1; |
---|
431 | tso.isNotMinimal = NULL; |
---|
432 | tso.p = syz; |
---|
433 | tso.order = pTotaldegree(syz); |
---|
434 | tso.syz = pHead(syz); |
---|
435 | pSetComp(tso.syz,newEl+1); |
---|
436 | tso.lcm = pHead(tso.syz); |
---|
437 | tso.length = pLength(syz); |
---|
438 | if (syOrder_Hilb(syz,syzstr,index,newEl+1)) |
---|
439 | syResetShiftedComponents(syzstr, index+1,1); |
---|
440 | rSetmS(tso.syz,syzstr->truecomponents[index],syzstr->ShiftedComponents[index]); |
---|
441 | #ifdef SHOW_PROT |
---|
442 | Print("erzeuge Halbpaar im Module %d,%d mit: \n",index,tso.order); |
---|
443 | Print("syz: ");pWrite(tso.syz); |
---|
444 | Print("sPoly: ");pWrite(tso.p); |
---|
445 | PrintLn(); |
---|
446 | #endif |
---|
447 | syEnterPair(syzstr,&tso,&l,index); |
---|
448 | } |
---|
449 | /*3 |
---|
450 | * computes the order of pairs of same degree |
---|
451 | * must be monoton |
---|
452 | */ |
---|
453 | static intvec* syLinStrat2(SSet nextPairs, syStrategy syzstr, |
---|
454 | int howmuch, int index,intvec ** secondpairs) |
---|
455 | { |
---|
456 | ideal o_r=syzstr->res[index+1]; |
---|
457 | int i=0,i1=0,i2=0,l,ll=IDELEMS(o_r); |
---|
458 | intvec *result=NewIntvec1(howmuch+1); |
---|
459 | BOOLEAN isDivisible; |
---|
460 | SObject tso; |
---|
461 | |
---|
462 | #ifndef USE_HEURISTIC2 |
---|
463 | while (i1<howmuch) |
---|
464 | { |
---|
465 | (*result)[i1] = i1+1; |
---|
466 | i1++; |
---|
467 | } |
---|
468 | return result; |
---|
469 | #else |
---|
470 | while ((ll>0) && (o_r->m[ll-1]==NULL)) ll--; |
---|
471 | while (i<howmuch) |
---|
472 | { |
---|
473 | tso = nextPairs[i]; |
---|
474 | isDivisible = FALSE; |
---|
475 | l = 0; |
---|
476 | while ((l<ll) && (!isDivisible)) |
---|
477 | { |
---|
478 | if (o_r->m[l]!=NULL) |
---|
479 | { |
---|
480 | isDivisible = isDivisible || |
---|
481 | pDivisibleBy1(o_r->m[l],tso.lcm); |
---|
482 | } |
---|
483 | l++; |
---|
484 | } |
---|
485 | if (isDivisible) |
---|
486 | { |
---|
487 | #ifdef SHOW_PROT |
---|
488 | Print("streiche Paar im Modul %d,%d mit: \n",index,nextPairs[i].order); |
---|
489 | Print("poly1: ");pWrite(nextPairs[i].p1); |
---|
490 | Print("poly2: ");pWrite(nextPairs[i].p2); |
---|
491 | Print("syz: ");pWrite(nextPairs[i].syz); |
---|
492 | Print("sPoly: ");pWrite(nextPairs[i].p); |
---|
493 | PrintLn(); |
---|
494 | #endif |
---|
495 | //syDeletePair(&nextPairs[i]); |
---|
496 | if (*secondpairs==NULL) *secondpairs = NewIntvec1(howmuch); |
---|
497 | (**secondpairs)[i2] = i+1; |
---|
498 | i2++; |
---|
499 | #ifdef SHOW_CRIT |
---|
500 | crit++; |
---|
501 | #endif |
---|
502 | } |
---|
503 | else |
---|
504 | { |
---|
505 | // nextPairs[i].p = sySPoly(tso.p1, tso.p2,tso.lcm); |
---|
506 | (*result)[i1] = i+1; |
---|
507 | i1++; |
---|
508 | } |
---|
509 | i++; |
---|
510 | } |
---|
511 | return result; |
---|
512 | #endif |
---|
513 | } |
---|
514 | |
---|
515 | inline void sySPRedSyz(syStrategy syzstr,sSObject redWith,poly q=NULL) |
---|
516 | { |
---|
517 | poly p=pDivide(q,redWith.p); |
---|
518 | pSetCoeff(p,nDiv(pGetCoeff(q),pGetCoeff(redWith.p))); |
---|
519 | int il=-1; |
---|
520 | kBucket_Minus_m_Mult_p(syzstr->syz_bucket,p,redWith.syz,&il,NULL); |
---|
521 | } |
---|
522 | |
---|
523 | static poly syRed_Hilb(poly toRed,syStrategy syzstr,int index) |
---|
524 | { |
---|
525 | ideal redWith=syzstr->res[index]; |
---|
526 | if (redWith==NULL) return toRed; |
---|
527 | int j=IDELEMS(redWith),i; |
---|
528 | poly q,result=NULL,resultp; |
---|
529 | |
---|
530 | while ((j>0) && (redWith->m[j-1]==NULL)) j--; |
---|
531 | if ((toRed==NULL) || (j==0)) return toRed; |
---|
532 | kBucketInit(syzstr->bucket,toRed,-1); |
---|
533 | q = kBucketGetLm(syzstr->bucket); |
---|
534 | loop |
---|
535 | { |
---|
536 | if (q==NULL) |
---|
537 | { |
---|
538 | break; |
---|
539 | } |
---|
540 | i = 0; |
---|
541 | loop |
---|
542 | { |
---|
543 | if (pDivisibleBy1(redWith->m[i],q)) |
---|
544 | { |
---|
545 | number up = kBucketPolyRed(syzstr->bucket,redWith->m[i], |
---|
546 | pLength(redWith->m[i]), NULL); |
---|
547 | nDelete(&up); |
---|
548 | q = kBucketGetLm(syzstr->bucket); |
---|
549 | if (toRed==NULL) break; |
---|
550 | i = 0; |
---|
551 | } |
---|
552 | else |
---|
553 | { |
---|
554 | i++; |
---|
555 | } |
---|
556 | if ((i>=j) || (q==NULL)) break; |
---|
557 | } |
---|
558 | if (q!=NULL) |
---|
559 | { |
---|
560 | if (result==NULL) |
---|
561 | { |
---|
562 | resultp = result = kBucketExtractLm(syzstr->bucket); |
---|
563 | } |
---|
564 | else |
---|
565 | { |
---|
566 | pNext(resultp) = kBucketExtractLm(syzstr->bucket); |
---|
567 | pIter(resultp); |
---|
568 | } |
---|
569 | q = kBucketGetLm(syzstr->bucket); |
---|
570 | } |
---|
571 | } |
---|
572 | kBucketClear(syzstr->bucket,&q,&i); |
---|
573 | if (q!=NULL) Print("Hier ist was schief gelaufen!\n"); |
---|
574 | return result; |
---|
575 | } |
---|
576 | |
---|
577 | intvec *ivStrip(intvec* arg) |
---|
578 | { |
---|
579 | int l=arg->rows()*arg->cols(),i=0,ii=0; |
---|
580 | intvec *tempV=NewIntvec1(l); |
---|
581 | |
---|
582 | while (i+ii<l) |
---|
583 | { |
---|
584 | if ((*arg)[i+ii]==0) |
---|
585 | { |
---|
586 | ii++; |
---|
587 | } |
---|
588 | else |
---|
589 | { |
---|
590 | (*tempV)[i] = (*arg)[i+ii]; |
---|
591 | i++; |
---|
592 | } |
---|
593 | } |
---|
594 | if (i==0) |
---|
595 | { |
---|
596 | delete tempV; |
---|
597 | return NULL; |
---|
598 | } |
---|
599 | intvec * result=NewIntvec1(i+1); |
---|
600 | for (ii=0;ii<i;ii++) |
---|
601 | (*result)[ii] = (*tempV)[ii]; |
---|
602 | delete tempV; |
---|
603 | return result; |
---|
604 | } |
---|
605 | |
---|
606 | /*3 |
---|
607 | * reduces all pairs of degree deg in the module index |
---|
608 | * put the reduced generators to the resolvente which contains |
---|
609 | * the truncated kStd |
---|
610 | */ |
---|
611 | static void syRedNextPairs_Hilb(SSet nextPairs, syStrategy syzstr, |
---|
612 | int howmuch, int index,int actord,int* toSub, |
---|
613 | int *maxindex,int *maxdeg) |
---|
614 | { |
---|
615 | int i,j,k=IDELEMS(syzstr->res[index]); |
---|
616 | int ks=IDELEMS(syzstr->res[index+1]),kk,l,ll; |
---|
617 | int ks1=IDELEMS(syzstr->orderedRes[index+1]); |
---|
618 | int kres=(*syzstr->Tl)[index]; |
---|
619 | int toGo=0; |
---|
620 | number coefgcd,n; |
---|
621 | SSet redset=syzstr->resPairs[index]; |
---|
622 | poly p=NULL,q,tp; |
---|
623 | intvec *spl1=NewIntvec1(howmuch+1); |
---|
624 | SObject tso; |
---|
625 | intvec *spl3=NULL; |
---|
626 | #ifdef USE_HEURISTIC1 |
---|
627 | intvec *spl2=NewIntvec3(howmuch+1,howmuch+1,0); |
---|
628 | int there_are_superfluous=0; |
---|
629 | int step=1,jj,j1,j2; |
---|
630 | #endif |
---|
631 | long * ShiftedComponents = syzstr->ShiftedComponents[index]; |
---|
632 | int* Components = syzstr->truecomponents[index]; |
---|
633 | assume(Components != NULL && ShiftedComponents != NULL); |
---|
634 | BOOLEAN need_reset; |
---|
635 | |
---|
636 | actord += index; |
---|
637 | if ((nextPairs==NULL) || (howmuch==0)) return; |
---|
638 | while ((k>0) && (syzstr->res[index]->m[k-1]==NULL)) k--; |
---|
639 | while ((ks>0) && (syzstr->res[index+1]->m[ks-1]==NULL)) ks--; |
---|
640 | while ((ks1>0) && (syzstr->orderedRes[index+1]->m[ks1-1]==NULL)) ks1--; |
---|
641 | while ((kres>0) && |
---|
642 | ((redset[kres-1].p==NULL) || (redset[kres-1].order>actord))) kres--; |
---|
643 | while ((kres<(*syzstr->Tl)[index]) && |
---|
644 | (redset[kres-1].order!=0) && (redset[kres-1].order<=actord)) kres++; |
---|
645 | spl1 = syLinStrat2(nextPairs,syzstr,howmuch,index,&spl3); |
---|
646 | #ifdef SHOW_PROT |
---|
647 | Print("spl1 ist hier: ");spl1->show(0,0); |
---|
648 | #endif |
---|
649 | i=0; |
---|
650 | kk = (*spl1)[i]-1; |
---|
651 | if (index==1) |
---|
652 | { |
---|
653 | intvec * temp1_hilb = hHstdSeries(syzstr->res[index],NULL,NULL,NULL); |
---|
654 | if (actord<temp1_hilb->length()) |
---|
655 | { |
---|
656 | toGo = (*temp1_hilb)[actord]; |
---|
657 | #ifdef SHOW_HILB |
---|
658 | Print("\nStze toGo im Modul %d und Grad %d auf: %d\n",1,actord-1,toGo); |
---|
659 | #endif |
---|
660 | } |
---|
661 | delete temp1_hilb; |
---|
662 | } |
---|
663 | else |
---|
664 | { |
---|
665 | if (actord<=(syzstr->hilb_coeffs[index])->length()) |
---|
666 | { |
---|
667 | toGo = (*syzstr->hilb_coeffs[index])[actord-1]; |
---|
668 | #ifdef SHOW_HILB |
---|
669 | Print("\nStze toGo im Modul %d und Grad %d auf: %d\n",index,actord-1,toGo); |
---|
670 | #endif |
---|
671 | } |
---|
672 | } |
---|
673 | if ((syzstr->hilb_coeffs[index+1]!=NULL) && |
---|
674 | (actord<=(syzstr->hilb_coeffs[index+1])->length())) |
---|
675 | { |
---|
676 | toGo += (*syzstr->hilb_coeffs[index+1])[actord-1]; |
---|
677 | #ifdef SHOW_HILB |
---|
678 | Print("\nAddiere zu toGo aus Modul %d und Grad %d: %d\n",index+1,actord-1,(*syzstr->hilb_coeffs[index+1])[actord-1]); |
---|
679 | #endif |
---|
680 | } |
---|
681 | #ifdef SHOW_HILB |
---|
682 | Print("<H%d>",toGo); |
---|
683 | #endif |
---|
684 | while (kk>=0) |
---|
685 | { |
---|
686 | if (toGo==0) |
---|
687 | { |
---|
688 | while (kk>=0) |
---|
689 | { |
---|
690 | pDelete(&nextPairs[kk].p); |
---|
691 | pDelete(&nextPairs[kk].syz); |
---|
692 | syDeletePair(&nextPairs[kk]); |
---|
693 | nextPairs[kk].p = nextPairs[kk].syz = nextPairs[kk].lcm = NULL; |
---|
694 | i++; |
---|
695 | kk = (*spl1)[i]-1; |
---|
696 | #ifdef USE_HEURISTIC2 |
---|
697 | if (kk<0) |
---|
698 | { |
---|
699 | i = 0; |
---|
700 | delete spl1; |
---|
701 | spl1 = spl3; |
---|
702 | spl3 = NULL; |
---|
703 | if (spl1!=NULL) |
---|
704 | kk = (*spl1)[i]-1; |
---|
705 | } |
---|
706 | #endif |
---|
707 | } |
---|
708 | if (spl1!=NULL) delete spl1; |
---|
709 | break; |
---|
710 | } |
---|
711 | tso = nextPairs[kk]; |
---|
712 | if ((tso.p1!=NULL) && (tso.p2!=NULL)) |
---|
713 | { |
---|
714 | #ifdef SHOW_CRIT |
---|
715 | cons_pairs++; |
---|
716 | #endif |
---|
717 | //tso.p = sySPoly(tso.p1, tso.p2,tso.lcm); |
---|
718 | tso.p = ksOldCreateSpoly(tso.p2, tso.p1); |
---|
719 | #ifdef SHOW_PROT |
---|
720 | Print("reduziere Paar mit: \n"); |
---|
721 | Print("poly1: ");pWrite(tso.p1); |
---|
722 | Print("poly2: ");pWrite(tso.p2); |
---|
723 | Print("syz: ");pWrite(tso.syz); |
---|
724 | Print("sPoly: ");pWrite(tso.p); |
---|
725 | #endif |
---|
726 | if (tso.p != NULL) |
---|
727 | { |
---|
728 | kBucketInit(syzstr->bucket,tso.p,-1); |
---|
729 | kBucketInit(syzstr->syz_bucket,tso.syz,-1); |
---|
730 | q = kBucketGetLm(syzstr->bucket); |
---|
731 | j = 0; |
---|
732 | while (j<kres) |
---|
733 | { |
---|
734 | if ((redset[j].p!=NULL) && (pDivisibleBy1(redset[j].p,q)) |
---|
735 | && ((redset[j].ind1!=tso.ind1) || (redset[j].ind2!=tso.ind2))) |
---|
736 | { |
---|
737 | #ifdef SHOW_RED |
---|
738 | Print("reduziere: ");pWrite(tso.p); |
---|
739 | Print("syz: ");pWrite(tso.syz); |
---|
740 | Print("mit: ");pWrite(redset[j].p); |
---|
741 | Print("syz: ");pWrite(redset[j].syz); |
---|
742 | #endif |
---|
743 | currcomponents = syzstr->truecomponents[index]; |
---|
744 | currShiftedComponents = syzstr->ShiftedComponents[index]; |
---|
745 | rChangeSComps(currcomponents, |
---|
746 | currShiftedComponents, |
---|
747 | IDELEMS(syzstr->res[index])); |
---|
748 | sySPRedSyz(syzstr,redset[j],q); |
---|
749 | currcomponents = syzstr->truecomponents[index-1]; |
---|
750 | currShiftedComponents = syzstr->ShiftedComponents[index-1]; |
---|
751 | rChangeSComps(currcomponents, |
---|
752 | currShiftedComponents, |
---|
753 | IDELEMS(syzstr->res[index-1])); |
---|
754 | number up = kBucketPolyRed(syzstr->bucket,redset[j].p, |
---|
755 | redset[j].length, NULL); |
---|
756 | nDelete(&up); |
---|
757 | q = kBucketGetLm(syzstr->bucket); |
---|
758 | #ifdef SHOW_RED |
---|
759 | Print("zu: ");pWrite(tso.p); |
---|
760 | Print("syz: ");pWrite(tso.syz); |
---|
761 | PrintLn(); |
---|
762 | #endif |
---|
763 | if (q==NULL) break; |
---|
764 | j = 0; |
---|
765 | } |
---|
766 | else |
---|
767 | { |
---|
768 | j++; |
---|
769 | } |
---|
770 | } |
---|
771 | kBucketClear(syzstr->bucket,&tso.p,&tso.length); |
---|
772 | currcomponents = syzstr->truecomponents[index]; |
---|
773 | currShiftedComponents = syzstr->ShiftedComponents[index]; |
---|
774 | rChangeSComps(currcomponents, |
---|
775 | currShiftedComponents, |
---|
776 | IDELEMS(syzstr->res[index])); |
---|
777 | int il; |
---|
778 | kBucketClear(syzstr->syz_bucket,&tso.syz,&il); |
---|
779 | currcomponents = syzstr->truecomponents[index-1]; |
---|
780 | currShiftedComponents = syzstr->ShiftedComponents[index-1]; |
---|
781 | rChangeSComps(currcomponents, |
---|
782 | currShiftedComponents, |
---|
783 | IDELEMS(syzstr->res[index-1])); |
---|
784 | } |
---|
785 | #ifdef SHOW_PROT |
---|
786 | Print("erhalte Paar mit: \n"); |
---|
787 | Print("syz: ");pWrite(tso.syz); |
---|
788 | Print("sPoly: ");pWrite(tso.p); |
---|
789 | PrintLn(); |
---|
790 | #endif |
---|
791 | #ifdef SHOW_SPRFL |
---|
792 | //PrintLn(); |
---|
793 | wrp(tso.lcm); |
---|
794 | Print(" mit index %d, %d ",tso.ind1,tso.ind2); |
---|
795 | #endif |
---|
796 | if (tso.p != NULL) |
---|
797 | { |
---|
798 | if (TEST_OPT_PROT) PrintS("g"); |
---|
799 | (*toSub)++; |
---|
800 | toGo--; |
---|
801 | if (!nIsOne(pGetCoeff(tso.p))) |
---|
802 | { |
---|
803 | number n=nInvers(pGetCoeff(tso.p)); |
---|
804 | pNorm(tso.p); |
---|
805 | pMultN(tso.syz,n); |
---|
806 | nDelete(&n); |
---|
807 | } |
---|
808 | if (k==IDELEMS((syzstr->res)[index])) |
---|
809 | syEnlargeFields(syzstr,index); |
---|
810 | syzstr->res[index]->m[k] = tso.p; |
---|
811 | k++; |
---|
812 | } |
---|
813 | else |
---|
814 | { |
---|
815 | if (ks==IDELEMS(syzstr->res[index+1])) |
---|
816 | syEnlargeFields(syzstr,index+1); |
---|
817 | currcomponents = syzstr->truecomponents[index]; |
---|
818 | currShiftedComponents = syzstr->ShiftedComponents[index]; |
---|
819 | rChangeSComps(currcomponents, |
---|
820 | currShiftedComponents, |
---|
821 | IDELEMS(syzstr->res[index])); |
---|
822 | syzstr->res[index+1]->m[ks] = syRed_Hilb(tso.syz,syzstr,index+1); |
---|
823 | currcomponents = syzstr->truecomponents[index-1]; |
---|
824 | currShiftedComponents = syzstr->ShiftedComponents[index-1]; |
---|
825 | rChangeSComps(currcomponents, |
---|
826 | currShiftedComponents, |
---|
827 | IDELEMS(syzstr->res[index-1])); |
---|
828 | if (syzstr->res[index+1]->m[ks]!=NULL) |
---|
829 | { |
---|
830 | if (TEST_OPT_PROT) PrintS("s"); |
---|
831 | toGo--; |
---|
832 | pNorm(syzstr->res[index+1]->m[ks]); |
---|
833 | syHalfPair(syzstr->res[index+1]->m[ks],ks1,syzstr,index+1); |
---|
834 | ks++; |
---|
835 | ks1++; |
---|
836 | if (index+1>*maxindex) *maxindex = index+1; |
---|
837 | if (actord-index>*maxdeg) *maxdeg = actord-index; |
---|
838 | } |
---|
839 | else |
---|
840 | { |
---|
841 | if (TEST_OPT_PROT) PrintS("-"); |
---|
842 | #ifdef SHOW_CRIT |
---|
843 | spfl++; |
---|
844 | #endif |
---|
845 | #ifdef USE_HEURISTIC1 |
---|
846 | if (there_are_superfluous>=0) |
---|
847 | { |
---|
848 | j = i+1; |
---|
849 | jj = (*spl1)[j]-1; |
---|
850 | j1 = 1; |
---|
851 | while (jj>=0) |
---|
852 | { |
---|
853 | if (tso.ind2==nextPairs[jj].ind2) |
---|
854 | { |
---|
855 | IMATELEM(*spl2,j1,step) = jj+1; |
---|
856 | j1++; |
---|
857 | for (j2=j;j2<spl1->length()-1;j2++) |
---|
858 | { |
---|
859 | (*spl1)[j2] = (*spl1)[j2+1]; |
---|
860 | } |
---|
861 | } |
---|
862 | else |
---|
863 | { |
---|
864 | j++; |
---|
865 | } |
---|
866 | jj = (*spl1)[j]-1; |
---|
867 | } |
---|
868 | step++; |
---|
869 | if (there_are_superfluous==0) there_are_superfluous = 1; |
---|
870 | } |
---|
871 | #endif |
---|
872 | #ifdef SHOW_SPRFL |
---|
873 | Print("ist ueberfluessig in Mod %d",index); |
---|
874 | //Print("\n ueberfluessig in Mod %d:",index); |
---|
875 | //wrp(tso.lcm); |
---|
876 | //PrintLn(); |
---|
877 | #endif |
---|
878 | } |
---|
879 | tso.syz = NULL; |
---|
880 | syDeletePair(&tso); |
---|
881 | tso.p = tso.syz = tso.lcm = NULL; |
---|
882 | } |
---|
883 | nextPairs[kk] = tso; |
---|
884 | } |
---|
885 | #ifdef SHOW_SPRFL |
---|
886 | PrintLn(); |
---|
887 | #endif |
---|
888 | i++; |
---|
889 | #ifdef SHOW_PROT |
---|
890 | Print("spl1 ist hier: ");spl1->show(0,0); |
---|
891 | Print("naechstes i ist: %d",i); |
---|
892 | #endif |
---|
893 | kk = (*spl1)[i]-1; |
---|
894 | #ifdef USE_HEURISTIC1 |
---|
895 | if ((kk<0) && (there_are_superfluous>0)) |
---|
896 | { |
---|
897 | i = 0; |
---|
898 | delete spl1; |
---|
899 | spl1 = ivStrip(spl2); |
---|
900 | delete spl2; |
---|
901 | if (spl1!=NULL) |
---|
902 | { |
---|
903 | there_are_superfluous = -1; |
---|
904 | kk = (*spl1)[i]-1; |
---|
905 | } |
---|
906 | } |
---|
907 | #endif |
---|
908 | #ifdef USE_HEURISTIC2 |
---|
909 | if ((kk<0) && (toGo>0)) |
---|
910 | { |
---|
911 | #ifdef SHOW_CRIT |
---|
912 | crit_fails++; |
---|
913 | #endif |
---|
914 | i = 0; |
---|
915 | delete spl1; |
---|
916 | spl1 = spl3; |
---|
917 | spl3 = NULL; |
---|
918 | if (spl1!=NULL) |
---|
919 | kk = (*spl1)[i]-1; |
---|
920 | } |
---|
921 | #endif |
---|
922 | } |
---|
923 | delete spl1; |
---|
924 | } |
---|
925 | |
---|
926 | void sySetNewHilb(syStrategy syzstr, int toSub,int index,int actord) |
---|
927 | { |
---|
928 | int i; |
---|
929 | actord += index; |
---|
930 | intvec * temp_hilb = hHstdSeries(syzstr->res[index+1],NULL,NULL,NULL); |
---|
931 | if (syzstr->hilb_coeffs[index+1]==NULL) |
---|
932 | { |
---|
933 | syzstr->hilb_coeffs[index+1] = NewIntvec1(16*((actord/16)+1)); |
---|
934 | } |
---|
935 | else if (actord>=syzstr->hilb_coeffs[index+1]->length()) |
---|
936 | { |
---|
937 | intvec * ttt=NewIntvec1(16*((actord/16)+1)); |
---|
938 | for (i=syzstr->hilb_coeffs[index+1]->length()-1;i>=0;i--) |
---|
939 | { |
---|
940 | (*ttt)[i] = (*(syzstr->hilb_coeffs[index+1]))[i]; |
---|
941 | } |
---|
942 | delete syzstr->hilb_coeffs[index+1]; |
---|
943 | syzstr->hilb_coeffs[index+1] = ttt; |
---|
944 | } |
---|
945 | if (actord+1<temp_hilb->length()) |
---|
946 | { |
---|
947 | #ifdef SHOW_HILB |
---|
948 | Print("\nSetze fuer Modul %d im Grad %d die Wert: \n",index+1,actord); |
---|
949 | (temp_hilb)->show(0,0); |
---|
950 | #endif |
---|
951 | int k=min(temp_hilb->length()-1,(syzstr->hilb_coeffs[index+1])->length()); |
---|
952 | for (int j=k;j>actord;j--) |
---|
953 | (*(syzstr->hilb_coeffs[index+1]))[j-1] = (*temp_hilb)[j]; |
---|
954 | } |
---|
955 | else |
---|
956 | { |
---|
957 | (*(syzstr->hilb_coeffs[index+1]))[actord] = 0; |
---|
958 | } |
---|
959 | delete temp_hilb; |
---|
960 | if ((index>1) && (actord<=syzstr->hilb_coeffs[index]->length())) |
---|
961 | { |
---|
962 | #ifdef SHOW_HILB |
---|
963 | Print("\nSubtrahiere im Modul %d im Grad %d den Wert: %d\n",index,actord-1,toSub); |
---|
964 | #endif |
---|
965 | (*syzstr->hilb_coeffs[index])[actord-1]-=toSub; |
---|
966 | } |
---|
967 | #ifdef SHOW_HILB |
---|
968 | Print("<h,%d>",(*(syzstr->hilb_coeffs[index+1]))[actord]); |
---|
969 | #endif |
---|
970 | } |
---|
971 | |
---|
972 | /*3 |
---|
973 | * reduces the generators of the module index in degree deg |
---|
974 | * (which are actual syzygies of the module index-1) |
---|
975 | * wrt. the ideal generated by elements of lower degrees |
---|
976 | */ |
---|
977 | static void syRedGenerOfCurrDeg_Hilb(syStrategy syzstr, int deg,int *maxindex,int *maxdeg) |
---|
978 | { |
---|
979 | ideal res=syzstr->res[1]; |
---|
980 | int i=0,j,k=IDELEMS(res),k1=IDELEMS(syzstr->orderedRes[1]); |
---|
981 | SSet sPairs1=syzstr->resPairs[1]; |
---|
982 | SSet sPairs=syzstr->resPairs[0]; |
---|
983 | |
---|
984 | while ((k>0) && (res->m[k-1]==NULL)) k--; |
---|
985 | while ((k1>0) && (syzstr->orderedRes[1]->m[k1-1]==NULL)) k1--; |
---|
986 | while ((i<(*syzstr->Tl)[0]) && (((sPairs)[i].syz==NULL) || |
---|
987 | ((sPairs)[i].order<deg))) |
---|
988 | i++; |
---|
989 | if ((i>=(*syzstr->Tl)[0]) || ((sPairs)[i].order>deg)) return; |
---|
990 | while ((i<(*syzstr->Tl)[0]) && (((sPairs)[i].syz==NULL) || |
---|
991 | ((sPairs)[i].order==deg))) |
---|
992 | { |
---|
993 | if ((sPairs)[i].syz!=NULL) |
---|
994 | { |
---|
995 | #ifdef SHOW_PROT |
---|
996 | Print("reduziere Erzeuger: \n"); |
---|
997 | Print("syz: ");pWrite((sPairs)[i].syz); |
---|
998 | #endif |
---|
999 | (sPairs)[i].syz = syRed_Hilb((sPairs)[i].syz,syzstr,1); |
---|
1000 | #ifdef SHOW_PROT |
---|
1001 | Print("erhalte Erzeuger: \n"); |
---|
1002 | Print("syz: ");pWrite((sPairs)[i].syz); |
---|
1003 | PrintLn(); |
---|
1004 | #endif |
---|
1005 | if ((sPairs)[i].syz != NULL) |
---|
1006 | { |
---|
1007 | if (k==IDELEMS(res)) |
---|
1008 | { |
---|
1009 | syEnlargeFields(syzstr,1); |
---|
1010 | res=syzstr->res[1]; |
---|
1011 | } |
---|
1012 | if (BTEST1(6)) |
---|
1013 | { |
---|
1014 | if ((sPairs)[i].isNotMinimal==NULL) |
---|
1015 | { |
---|
1016 | PrintS("\nminimal generator: "); |
---|
1017 | pWrite((syzstr->resPairs[0])[i].syz); |
---|
1018 | PrintS("comes from: ");pWrite((syzstr->resPairs[0])[i].p1); |
---|
1019 | PrintS("and: ");pWrite((syzstr->resPairs[0])[i].p2); |
---|
1020 | } |
---|
1021 | } |
---|
1022 | res->m[k] = (sPairs)[i].syz; |
---|
1023 | pNorm(res->m[k]); |
---|
1024 | syHalfPair(res->m[k],k1,syzstr,1); |
---|
1025 | k1++; |
---|
1026 | k++; |
---|
1027 | if (1>*maxindex) *maxindex = 1; |
---|
1028 | if (deg-1>*maxdeg) *maxdeg = deg-1; |
---|
1029 | } |
---|
1030 | } |
---|
1031 | i++; |
---|
1032 | } |
---|
1033 | } |
---|
1034 | |
---|
1035 | /*3 |
---|
1036 | * reorders the result (stored in orderedRes) according |
---|
1037 | * to the seqence given by res |
---|
1038 | */ |
---|
1039 | static void syReOrdResult_Hilb(syStrategy syzstr,int maxindex,int maxdeg) |
---|
1040 | { |
---|
1041 | ideal reor,toreor; |
---|
1042 | int i,j,k,l,m,togo; |
---|
1043 | syzstr->betti = NewIntvec3(maxdeg,maxindex+1,0); |
---|
1044 | (*syzstr->betti)[0] = 1; |
---|
1045 | for (i=1;i<=syzstr->length;i++) |
---|
1046 | { |
---|
1047 | if (!idIs0(syzstr->orderedRes[i])) |
---|
1048 | { |
---|
1049 | toreor = syzstr->orderedRes[i]; |
---|
1050 | k = IDELEMS(toreor); |
---|
1051 | while ((k>0) && (toreor->m[k-1]==NULL)) k--; |
---|
1052 | reor = idInit(k,toreor->rank); |
---|
1053 | togo = IDELEMS(syzstr->res[i]); |
---|
1054 | for (j=0;j<k;j++) |
---|
1055 | { |
---|
1056 | if (toreor->m[j]!=NULL) (IMATELEM(*syzstr->betti,pFDeg(toreor->m[j])-i+1,i+1))++; |
---|
1057 | reor->m[j] = toreor->m[j]; |
---|
1058 | toreor->m[j] = NULL; |
---|
1059 | } |
---|
1060 | m = 0; |
---|
1061 | for (j=0;j<togo;j++) |
---|
1062 | { |
---|
1063 | if (syzstr->res[i]->m[j]!=NULL) |
---|
1064 | { |
---|
1065 | l = 0; |
---|
1066 | while ((l<k) && (syzstr->res[i]->m[j]!=reor->m[l])) l++; |
---|
1067 | if (l<k) |
---|
1068 | { |
---|
1069 | toreor->m[m] = reor->m[l]; |
---|
1070 | reor->m[l] = NULL; |
---|
1071 | m++; |
---|
1072 | } |
---|
1073 | } |
---|
1074 | } |
---|
1075 | idDelete(&reor); |
---|
1076 | } |
---|
1077 | } |
---|
1078 | } |
---|
1079 | |
---|
1080 | /*2 |
---|
1081 | * the CoCoA-algorithm for free resolutions, using a formula |
---|
1082 | * for remaining pairs based on Hilbert-functions |
---|
1083 | */ |
---|
1084 | syStrategy syHilb(ideal arg,int * length) |
---|
1085 | { |
---|
1086 | int i,j,*ord,*b0,*b1,actdeg=32000,index=0,reg=-1; |
---|
1087 | int startdeg,howmuch,toSub=0; |
---|
1088 | int maxindex=0,maxdeg=0; |
---|
1089 | ideal temp=NULL; |
---|
1090 | SSet nextPairs; |
---|
1091 | ring tmpR=NULL; |
---|
1092 | ring origR = currRing; |
---|
1093 | pSetmProc oldSetm=pSetm; |
---|
1094 | pCompProc oldComp0=pComp0; |
---|
1095 | syStrategy syzstr=(syStrategy)Alloc0SizeOf(ssyStrategy); |
---|
1096 | |
---|
1097 | if ((idIs0(arg)) || (idRankFreeModule(arg)>0)) |
---|
1098 | { |
---|
1099 | syzstr->minres = (resolvente)Alloc0(sizeof(ideal)); |
---|
1100 | syzstr->length = 1; |
---|
1101 | syzstr->minres[0] = idInit(1,arg->rank); |
---|
1102 | return syzstr; |
---|
1103 | } |
---|
1104 | tmpR = (ring) Alloc0SizeOf(sip_sring); |
---|
1105 | tmpR->wvhdl = (int **)Alloc0(3 * sizeof(int *)); |
---|
1106 | ord = (int*)Alloc0(3*sizeof(int)); |
---|
1107 | b0 = (int*)Alloc0(3*sizeof(int)); |
---|
1108 | b1 = (int*)Alloc0(3*sizeof(int)); |
---|
1109 | ord[0] = ringorder_dp; |
---|
1110 | ord[1] = ringorder_S; |
---|
1111 | b0[0] = 1; |
---|
1112 | b1[0] = currRing->N; |
---|
1113 | tmpR->OrdSgn = 1; |
---|
1114 | tmpR->N = currRing->N; |
---|
1115 | tmpR->ch = currRing->ch; |
---|
1116 | tmpR->order = ord; |
---|
1117 | tmpR->block0 = b0; |
---|
1118 | tmpR->block1 = b1; |
---|
1119 | tmpR->P = currRing->P; |
---|
1120 | if (currRing->parameter!=NULL) |
---|
1121 | { |
---|
1122 | tmpR->minpoly=nCopy(currRing->minpoly); |
---|
1123 | tmpR->parameter=(char **)Alloc(rPar(currRing)*sizeof(char *)); |
---|
1124 | for(i=0;i<tmpR->P;i++) |
---|
1125 | { |
---|
1126 | tmpR->parameter[i]=mstrdup(currRing->parameter[i]); |
---|
1127 | } |
---|
1128 | } |
---|
1129 | tmpR->names = (char **)Alloc(tmpR->N * sizeof(char *)); |
---|
1130 | for (i=0; i<tmpR->N; i++) |
---|
1131 | { |
---|
1132 | tmpR->names[i] = mstrdup(currRing->names[i]); |
---|
1133 | } |
---|
1134 | currcomponents = (int*)Alloc0((arg->rank+1)*sizeof(int)); |
---|
1135 | currShiftedComponents = (long*)Alloc0((arg->rank+1)*sizeof(long)); |
---|
1136 | for (i=0;i<=arg->rank;i++) |
---|
1137 | { |
---|
1138 | currShiftedComponents[i] = (i)*SYZ_SHIFT_BASE; |
---|
1139 | currcomponents[i] = i; |
---|
1140 | } |
---|
1141 | rComplete(tmpR, 1); |
---|
1142 | rChangeCurrRing(tmpR, TRUE); |
---|
1143 | rChangeSComps(currcomponents, currShiftedComponents, arg->rank); |
---|
1144 | syzstr->syRing = tmpR; |
---|
1145 | /*--- initializes the data structures---------------*/ |
---|
1146 | #ifdef SHOW_CRIT |
---|
1147 | crit = 0; |
---|
1148 | crit1 = 0; |
---|
1149 | spfl = 0; |
---|
1150 | cons_pairs = 0; |
---|
1151 | crit_fails = 0; |
---|
1152 | #endif |
---|
1153 | syzstr->length = *length = pVariables+2; |
---|
1154 | syzstr->Tl = NewIntvec1(*length+1); |
---|
1155 | temp = idInit(IDELEMS(arg),arg->rank); |
---|
1156 | for (i=0;i<IDELEMS(arg);i++) |
---|
1157 | { |
---|
1158 | temp->m[i] = pFetchCopy(origR, arg->m[i]); |
---|
1159 | if (temp->m[i]!=NULL) |
---|
1160 | { |
---|
1161 | j = pTotaldegree(temp->m[i]); |
---|
1162 | if (j<actdeg) actdeg = j; |
---|
1163 | } |
---|
1164 | } |
---|
1165 | idTest(temp); |
---|
1166 | idSkipZeroes(temp); |
---|
1167 | syzstr->resPairs = syInitRes(temp,length,syzstr->Tl,syzstr->cw); |
---|
1168 | Free((ADDRESS)currcomponents,(arg->rank+1)*sizeof(int)); |
---|
1169 | Free((ADDRESS)currShiftedComponents,(arg->rank+1)*sizeof(int)); |
---|
1170 | syzstr->res = (resolvente)Alloc0((*length+1)*sizeof(ideal)); |
---|
1171 | syzstr->orderedRes = (resolvente)Alloc0((*length+1)*sizeof(ideal)); |
---|
1172 | syzstr->elemLength = (int**)Alloc0((*length+1)*sizeof(int*)); |
---|
1173 | syzstr->truecomponents = (int**)Alloc0((*length+1)*sizeof(int*)); |
---|
1174 | syzstr->backcomponents = (int**)Alloc0((*length+1)*sizeof(int*)); |
---|
1175 | syzstr->ShiftedComponents = (long**)Alloc0((*length+1)*sizeof(long*)); |
---|
1176 | syzstr->Howmuch = (int**)Alloc0((*length+1)*sizeof(int*)); |
---|
1177 | syzstr->Firstelem = (int**)Alloc0((*length+1)*sizeof(int*)); |
---|
1178 | syzstr->hilb_coeffs = (intvec**)Alloc0((*length+1)*sizeof(intvec*)); |
---|
1179 | syzstr->sev = (unsigned long **)Alloc0((*length+1)*sizeof(unsigned long*)); |
---|
1180 | syzstr->bucket = kBucketCreate(); |
---|
1181 | syzstr->syz_bucket = kBucketCreate(); |
---|
1182 | startdeg = actdeg; |
---|
1183 | nextPairs = syChosePairs(syzstr,&index,&howmuch,&actdeg); |
---|
1184 | /*--- computes the resolution ----------------------*/ |
---|
1185 | while (nextPairs!=NULL) |
---|
1186 | { |
---|
1187 | #ifdef SHOW_PROT |
---|
1188 | Print("compute %d Paare im Module %d im Grad %d \n",howmuch,index,actdeg+index); |
---|
1189 | #endif |
---|
1190 | if (TEST_OPT_PROT) Print("%d",actdeg); |
---|
1191 | if (TEST_OPT_PROT) Print("(m%d)",index); |
---|
1192 | if (index==0) |
---|
1193 | i = syInitSyzMod(syzstr,index,idRankFreeModule(arg)+1); |
---|
1194 | else |
---|
1195 | i = syInitSyzMod(syzstr,index); |
---|
1196 | currcomponents = syzstr->truecomponents[max(index-1,0)]; |
---|
1197 | currShiftedComponents = syzstr->ShiftedComponents[max(index-1,0)]; |
---|
1198 | rChangeSComps(currcomponents, currShiftedComponents, |
---|
1199 | IDELEMS(syzstr->res[max(index-1,0)])); |
---|
1200 | j = syInitSyzMod(syzstr,index+1); |
---|
1201 | if (index>0) |
---|
1202 | { |
---|
1203 | syRedNextPairs_Hilb(nextPairs,syzstr,howmuch,index,actdeg,&toSub,&maxindex,&maxdeg); |
---|
1204 | sySetNewHilb(syzstr,toSub,index,actdeg); |
---|
1205 | toSub = 0; |
---|
1206 | syCompactifyPairSet(syzstr->resPairs[index],(*syzstr->Tl)[index],0); |
---|
1207 | } |
---|
1208 | else |
---|
1209 | syRedGenerOfCurrDeg_Hilb(syzstr,actdeg,&maxindex,&maxdeg); |
---|
1210 | /*--- creates new pairs -----------------------------*/ |
---|
1211 | #ifdef SHOW_PROT |
---|
1212 | Print("Bilde neue Paare in Modul %d!\n",index); |
---|
1213 | #endif |
---|
1214 | syCreateNewPairs_Hilb(syzstr,index,actdeg); |
---|
1215 | if (index<(*length)-1) |
---|
1216 | { |
---|
1217 | #ifdef SHOW_PROT |
---|
1218 | Print("Bilde neue Paare in Modul %d!\n",index+1); |
---|
1219 | #endif |
---|
1220 | syCreateNewPairs_Hilb(syzstr,index+1,actdeg-1); } |
---|
1221 | index++; |
---|
1222 | nextPairs = syChosePairs(syzstr,&index,&howmuch,&actdeg); |
---|
1223 | } |
---|
1224 | syReOrdResult_Hilb(syzstr,maxindex,maxdeg); |
---|
1225 | #ifdef SHOW_RESULT |
---|
1226 | Print("minimal resolution:\n"); |
---|
1227 | for (int ti=1;ti<=*length;ti++) |
---|
1228 | { |
---|
1229 | if (!idIs0(syzstr->orderedRes[ti])) idPrint(syzstr->orderedRes[ti]); |
---|
1230 | } |
---|
1231 | Print("full resolution:\n"); |
---|
1232 | for (int ti=1;ti<=*length;ti++) |
---|
1233 | { |
---|
1234 | if (!idIs0(syzstr->res[ti])) idPrint(syzstr->res[ti]); |
---|
1235 | } |
---|
1236 | #endif |
---|
1237 | #ifdef SHOW_CRIT |
---|
1238 | Print("Criterion %d times applied\n",crit); |
---|
1239 | Print("Criterion1 %d times applied\n",crit1); |
---|
1240 | Print("%d superfluous pairs\n",spfl); |
---|
1241 | Print("%d pairs considered\n",cons_pairs); |
---|
1242 | Print("Criterion fails %d times\n",crit_fails); |
---|
1243 | crit = 0; |
---|
1244 | crit1 = 0; |
---|
1245 | spfl = 0; |
---|
1246 | cons_pairs = 0; |
---|
1247 | crit_fails = 0; |
---|
1248 | #endif |
---|
1249 | if (temp!=NULL) idDelete(&temp); |
---|
1250 | kBucketDestroy(&(syzstr->bucket)); |
---|
1251 | kBucketDestroy(&(syzstr->syz_bucket)); |
---|
1252 | if (ord!=NULL) |
---|
1253 | { |
---|
1254 | rChangeCurrRing(origR,TRUE); |
---|
1255 | } |
---|
1256 | else |
---|
1257 | { |
---|
1258 | pSetm=oldSetm; |
---|
1259 | pComp0=oldComp0; |
---|
1260 | } |
---|
1261 | if (TEST_OPT_PROT) PrintLn(); |
---|
1262 | return syzstr; |
---|
1263 | } |
---|