My Project
Loading...
Searching...
No Matches
MinorProcessor.cc
Go to the documentation of this file.
1
2
3
4#include "kernel/mod2.h"
5
7
8#include "polys/kbuckets.h"
9
10#include "kernel/structs.h"
11#include "kernel/polys.h"
13
14#include "kernel/ideals.h"
15
16using namespace std;
17
18#ifdef COUNT_AND_PRINT_OPERATIONS
19VAR long addsPoly = 0; /* for the number of additions of two polynomials */
20VAR long multsPoly = 0; /* for the number of multiplications of two polynomials */
21VAR long addsPolyForDiv = 0; /* for the number of additions of two polynomials for
22 polynomial division part */
23VAR long multsPolyForDiv = 0; /* for the number of multiplications of two polynomials
24 for polynomial division part */
25VAR long multsMon = 0; /* for the number of multiplications of two monomials */
26VAR long multsMonForDiv = 0; /* for the number of m-m-multiplications for polynomial
27 division part */
28VAR long savedMultsMFD = 0; /* number of m-m-multiplications that could be saved
29 when polynomial division would be optimal
30 (if p / t1 = t2 + ..., then t1 * t2 = LT(p), i.e.,
31 this multiplication need not be performed which
32 would save one m-m-multiplication) */
33VAR long divsMon = 0; /* for the number of divisions of two monomials;
34 these are all guaranteed to work, i.e., m1/m2 only
35 when exponentVector(m1) >= exponentVector(m2) */
36void printCounters (char* prefix, bool resetToZero)
37{
38 printf("%s [p+p(div) | p*p(div) | m*m(div, -save) | m/m ]", prefix);
39 printf(" = [%ld(%ld) | %ld(%ld) | %ld(%d, -%ld) | %ld]\n",
42 if (resetToZero)
43 {
44 multsMon = 0; addsPoly = 0; multsPoly = 0; divsMon = 0;
47 }
48}
49#endif
50/* COUNT_AND_PRINT_OPERATIONS */
51
53{
54 PrintS(this->toString().c_str());
55}
56
57int MinorProcessor::getBestLine (const int k, const MinorKey& mk) const
58{
59 /* This method identifies the row or column with the most zeros.
60 The returned index (bestIndex) is absolute within the pre-
61 defined matrix.
62 If some row has the most zeros, then the absolute (0-based)
63 row index is returned.
64 If, contrariwise, some column has the most zeros, then -1 minus
65 the absolute (0-based) column index is returned. */
66 int numberOfZeros = 0;
67 int bestIndex = 100000; /* We start with an invalid row/column index. */
68 int maxNumberOfZeros = -1; /* We update this variable whenever we find
69 a new so-far optimal row or column. */
70 for (int r = 0; r < k; r++)
71 {
72 /* iterate through all k rows of the momentary minor */
73 int absoluteR = mk.getAbsoluteRowIndex(r);
74 numberOfZeros = 0;
75 for (int c = 0; c < k; c++)
76 {
77 int absoluteC = mk.getAbsoluteColumnIndex(c);
79 }
81 {
82 /* We found a new best line which is a row. */
85 }
86 };
87 for (int c = 0; c < k; c++)
88 {
89 int absoluteC = mk.getAbsoluteColumnIndex(c);
90 numberOfZeros = 0;
91 for (int r = 0; r < k; r++)
92 {
93 int absoluteR = mk.getAbsoluteRowIndex(r);
95 }
97 {
98 /* We found a new best line which is a column. So we transform
99 the return value. Note that we can easily retrieve absoluteC
100 from bestLine: absoluteC = - 1 - bestLine. */
101 bestIndex = - absoluteC - 1;
103 }
104 };
105 return bestIndex;
106}
107
113
118
123
128
130 const int* rowIndices,
131 const int numberOfColumns,
132 const int* columnIndices)
133{
134 /* The method assumes ascending row and column indices in the
135 two argument arrays. These indices are understood to be zero-based.
136 The method will set the two arrays of ints in _container.
137 Example: The indices 0, 2, 3, 7 will be converted to an array with
138 one int representing the binary number 10001101
139 (check bits from right to left). */
140
143 int rowBlockCount = (highestRowIndex / 32) + 1;
144 unsigned *rowBlocks=(unsigned*)omAlloc(rowBlockCount*sizeof(unsigned));
145 for (int i = 0; i < rowBlockCount; i++) rowBlocks[i] = 0;
146 for (int i = 0; i < numberOfRows; i++)
147 {
148 int blockIndex = rowIndices[i] / 32;
149 int offset = rowIndices[i] % 32;
150 rowBlocks[blockIndex] += (1 << offset);
151 }
152
155 int columnBlockCount = (highestColumnIndex / 32) + 1;
156 unsigned *columnBlocks=(unsigned*)omAlloc0(columnBlockCount*sizeof(unsigned));
157 for (int i = 0; i < numberOfColumns; i++)
158 {
159 int blockIndex = columnIndices[i] / 32;
160 int offset = columnIndices[i] % 32;
161 columnBlocks[blockIndex] += (1 << offset);
162 }
163
167}
168
170{
171 /* This method moves _minor to the next valid (k x k)-minor within
172 _container. It returns true iff this is successful, i.e. iff
173 _minor did not already encode the terminal (k x k)-minor. */
174 if (_minor.compare(MinorKey(0, 0, 0, 0)) == 0)
175 {
176 /* This means that we haven't started yet. Thus, we are about
177 to compute the first (k x k)-minor. */
180 return true;
181 }
183 {
184 /* Here we were able to pick a next subset of columns
185 within the same subset of rows. */
186 return true;
187 }
189 {
190 /* Here we were not able to pick a next subset of columns
191 within the same subset of rows. But we could pick a next
192 subset of rows. We must hence reset the subset of columns: */
194 return true;
195 }
196 else
197 {
198 /* We were neither able to pick a next subset
199 of columns nor of rows. I.e., we have iterated through
200 all sensible choices of subsets of rows and columns. */
201 return false;
202 }
203}
204
205bool MinorProcessor::isEntryZero (const int /*absoluteRowIndex*/,
206 const int /*absoluteColumnIndex*/) const
207{
208 assume(false);
209 return false;
210}
211
213{
214 assume(false);
215 return "";
216}
217
218int MinorProcessor::IOverJ(const int i, const int j)
219{
220 /* This is a non-recursive implementation. */
221 assume( (i >= 0) && (j >= 0) && (i >= j));
222 if (j == 0 || i == j) return 1;
223 #if 0
224 int result = 1;
225 for (int k = i - j + 1; k <= i; k++) result *= k;
226 /* Now, result = (i - j + 1) * ... * i. */
227 for (int k = 2; k <= j; k++) result /= k;
228 /* Now, result = (i - j + 1) * ... * i / 1 / 2 ...
229 ... / j = i! / j! / (i - j)!. */
230 return result;
231 #endif
232 return binom(i,j);
233}
234
236{
237 /* This is a non-recursive implementation. */
238 assume(i >= 0);
239 int result = 1;
240 for (int j = 1; j <= i; j++) result *= j;
241 // Now, result = 1 * 2 * ... * i = i!
242 return result;
243}
244
245int MinorProcessor::NumberOfRetrievals (const int rows, const int columns,
246 const int containerMinorSize,
247 const int minorSize,
248 const bool multipleMinors)
249{
250 /* This method computes the number of potential retrievals
251 of a single minor when computing all minors of a given size
252 within a matrix of given size. */
253 int result = 0;
254 if (multipleMinors)
255 {
256 /* Here, we would like to compute all minors of size
257 containerMinorSize x containerMinorSize in a matrix
258 of size rows x columns.
259 Then, we need to retrieve any minor of size
260 minorSize x minorSize exactly n times, where n is as
261 follows: */
265 }
266 else
267 {
268 /* Here, we would like to compute just one minor of size
269 containerMinorSize x containerMinorSize. Then, we need
270 to retrieve any minor of size minorSize x minorSize exactly
271 (containerMinorSize - minorSize)! times: */
273 }
274 return result;
275}
276
278 _container(0, NULL, 0, NULL),
279 _minor(0, NULL, 0, NULL),
280 _containerRows(0),
281 _containerColumns(0),
282 _minorSize(0),
283 _rows(0),
284 _columns(0)
285{
286}
287
289
294
296{
297 char h[32];
298 string t = "";
299 string s = "IntMinorProcessor:";
300 s += "\n matrix: ";
301 sprintf(h, "%d", _rows); s += h;
302 s += " x ";
303 sprintf(h, "%d", _columns); s += h;
304 for (int r = 0; r < _rows; r++)
305 {
306 s += "\n ";
307 for (int c = 0; c < _columns; c++)
308 {
309 sprintf(h, "%d", getEntry(r, c)); t = h;
310 for (int k = 0; k < int(4 - strlen(h)); k++) s += " ";
311 s += t;
312 }
313 }
314 int myIndexArray[500];
315 s += "\n considered submatrix has row indices: ";
317 for (int k = 0; k < _containerRows; k++)
318 {
319 if (k != 0) s += ", ";
320 sprintf(h, "%d", myIndexArray[k]); s += h;
321 }
322 s += " (first row of matrix has index 0)";
323 s += "\n considered submatrix has column indices: ";
325 for (int k = 0; k < _containerColumns; k++)
326 {
327 if (k != 0) s += ", ";
328 sprintf(h, "%d", myIndexArray[k]); s += h;
329 }
330 s += " (first column of matrix has index 0)";
331 s += "\n size of considered minor(s): ";
332 sprintf(h, "%d", _minorSize); s += h;
333 s += "x";
334 s += h;
335 return s;
336}
337
339{
340 /* free memory of _intMatrix */
341 delete [] _intMatrix; _intMatrix = 0;
342}
343
349
351 const int numberOfColumns,
352 const int* matrix)
353{
354 /* free memory of _intMatrix */
356
359
360 /* allocate memory for new entries in _intMatrix */
361 int n = _rows * _columns;
362 _intMatrix = (int*)omAlloc(n*sizeof(int));
363
364 /* copying values from one-dimensional method
365 parameter "matrix" */
366 for (int i = 0; i < n; i++)
367 _intMatrix[i] = matrix[i];
368}
369
371 const int* rowIndices,
372 const int* columnIndices,
374 const int characteristic,
375 const ideal& iSB)
376{
379 /* call a helper method which recursively computes the minor using the
380 cache c: */
383}
384
386 const int* rowIndices,
387 const int* columnIndices,
388 const int characteristic,
389 const ideal& iSB,
390 const char* algorithm)
391{
394
395 /* call a helper method which computes the minor (without a cache): */
396 if (strcmp(algorithm, "Laplace") == 0)
398 iSB);
399 else if (strcmp(algorithm, "Bareiss") == 0)
401 iSB);
402 else assume(false);
403
404 /* The following code is never reached and just there to make the
405 compiler happy: */
406 return IntMinorValue();
407}
408
410 const ideal& iSB,
411 const char* algorithm)
412{
413 /* call a helper method which computes the minor (without a cache): */
414 if (strcmp(algorithm, "Laplace") == 0)
416 else if (strcmp(algorithm, "Bareiss") == 0)
418 else assume(false);
419
420 /* The following code is never reached and just there to make the
421 compiler happy: */
422 return IntMinorValue();
423}
424
426 IntMinorValue>& c,
427 const int characteristic,
428 const ideal& iSB)
429{
430 /* computation with cache */
432 iSB);
433}
434
435/* computes the reduction of an integer i modulo an ideal
436 which captures a std basis */
437int getReduction (const int i, const ideal& iSB)
438{
439 if (i == 0) return 0;
440 poly f = pISet(i);
441 poly g = kNF(iSB, currRing->qideal, f);
442 int result = 0;
443 if (g != NULL) result = n_Int(pGetCoeff(g), currRing->cf);
444 pDelete(&f);
445 pDelete(&g);
446 return result;
447}
448
450 const int k,
451 const MinorKey& mk,
452 const int characteristic,
453 const ideal& iSB)
454{
455 assume(k > 0); /* k is the minor's dimension; the minor must be at least
456 1x1 */
457 /* The method works by recursion, and using Lapace's Theorem along the
458 row/column with the most zeros. */
459 if (k == 1)
460 {
461 int e = getEntry(mk.getAbsoluteRowIndex(0), mk.getAbsoluteColumnIndex(0));
462 if (characteristic != 0) e = e % characteristic;
463 if (iSB != 0) e = getReduction(e, iSB);
464 return IntMinorValue(e, 0, 0, 0, 0, -1, -1); /* "-1" is to signal that any
465 statistics about the number
466 of retrievals does not make
467 sense, as we do not use a
468 cache. */
469 }
470 else
471 {
472 /* Here, the minor must be 2x2 or larger. */
473 int b = getBestLine(k, mk); /* row or column with most
474 zeros */
475 int result = 0; /* This will contain the
476 value of the minor. */
477 int s = 0; int m = 0; int as = 0; int am = 0; /* counters for additions and
478 multiplications, ..."a*"
479 for accumulated operation
480 counters */
481 bool hadNonZeroEntry = false;
482 if (b >= 0)
483 {
484 /* This means that the best line is the row with absolute (0-based)
485 index b.
486 Using Laplace, the sign of the contributing minors must be iterating;
487 the initial sign depends on the relative index of b in minorRowKey: */
488 int sign = (mk.getRelativeRowIndex(b) % 2 == 0 ? 1 : -1);
489 for (int c = 0; c < k; c++) /* This iterates over all involved columns. */
490 {
491 int absoluteC = mk.getAbsoluteColumnIndex(c);
492 if (getEntry(b, absoluteC) != 0) /* Only then do we have to consider
493 this sub-determinante. */
494 {
495 hadNonZeroEntry = true;
496 /* Next MinorKey is mk with row b and column absoluteC omitted: */
497 MinorKey subMk = mk.getSubMinorKey(b, absoluteC);
499 characteristic, iSB); /* recursive call */
500 m += mv.getMultiplications();
501 s += mv.getAdditions();
502 am += mv.getAccumulatedMultiplications();
503 as += mv.getAccumulatedAdditions();
504 /* adding sub-determinante times matrix entry
505 times appropriate sign: */
506 result += sign * mv.getResult() * getEntry(b, absoluteC);
507
509 s++; m++; as++, am++; /* This is for the last addition and
510 multiplication. */
511 }
512 sign = - sign; /* alternating the sign */
513 }
514 }
515 else
516 {
517 b = - b - 1;
518 /* This means that the best line is the column with absolute (0-based)
519 index b.
520 Using Laplace, the sign of the contributing minors must be iterating;
521 the initial sign depends on the relative index of b in
522 minorColumnKey: */
523 int sign = (mk.getRelativeColumnIndex(b) % 2 == 0 ? 1 : -1);
524 for (int r = 0; r < k; r++) /* This iterates over all involved rows. */
525 {
526 int absoluteR = mk.getAbsoluteRowIndex(r);
527 if (getEntry(absoluteR, b) != 0) /* Only then do we have to consider
528 this sub-determinante. */
529 {
530 hadNonZeroEntry = true;
531 /* Next MinorKey is mk with row absoluteR and column b omitted. */
532 MinorKey subMk = mk.getSubMinorKey(absoluteR, b);
533 IntMinorValue mv = getMinorPrivateLaplace(k - 1, subMk, characteristic, iSB); /* recursive call */
534 m += mv.getMultiplications();
535 s += mv.getAdditions();
536 am += mv.getAccumulatedMultiplications();
537 as += mv.getAccumulatedAdditions();
538 /* adding sub-determinante times matrix entry
539 times appropriate sign: */
540 result += sign * mv.getResult() * getEntry(absoluteR, b);
542 s++; m++; as++, am++; /* This is for the last addition and
543 multiplication. */
544 }
545 sign = - sign; /* alternating the sign */
546 }
547 }
548 if (hadNonZeroEntry)
549 {
550 s--; as--; /* first addition was 0 + ..., so we do not count it */
551 }
552 if (s < 0) s = 0; /* may happen when all subminors are zero and no
553 addition needs to be performed */
554 if (as < 0) as = 0; /* may happen when all subminors are zero and no
555 addition needs to be performed */
556 if (iSB != 0) result = getReduction(result, iSB);
557 IntMinorValue newMV(result, m, s, am, as, -1, -1);
558 /* "-1" is to signal that any statistics about the number of retrievals
559 does not make sense, as we do not use a cache. */
560 return newMV;
561 }
562}
563
564/* This method can only be used in the case of coefficients
565 coming from a field or at least from an integral domain. */
567 const int k,
568 const MinorKey& mk,
569 const int characteristic,
570 const ideal& iSB)
571{
572 assume(k > 0); /* k is the minor's dimension; the minor must be at least
573 1x1 */
574 int *theRows=(int*)omAlloc(k*sizeof(int));
575 mk.getAbsoluteRowIndices(theRows);
576 int *theColumns=(int*)omAlloc(k*sizeof(int));
577 mk.getAbsoluteColumnIndices(theColumns);
578 /* the next line provides the return value for the case k = 1 */
579 int e = getEntry(theRows[0], theColumns[0]);
580 if (characteristic != 0) e = e % characteristic;
581 if (iSB != 0) e = getReduction(e, iSB);
582 IntMinorValue mv(e, 0, 0, 0, 0, -1, -1);
583 if (k > 1)
584 {
585 /* the matrix to perform Bareiss with */
586 long *tempMatrix=(long*)omAlloc(k * k *sizeof(long));
587 /* copy correct set of entries from _intMatrix to tempMatrix */
588 int i = 0;
589 for (int r = 0; r < k; r++)
590 for (int c = 0; c < k; c++)
591 {
592 e = getEntry(theRows[r], theColumns[c]);
593 if (characteristic != 0) e = e % characteristic;
594 tempMatrix[i++] = e;
595 }
596 /* Bareiss algorithm operating on tempMatrix which is at least 2x2 */
597 int sign = 1; /* This will store the correct sign resulting
598 from permuting the rows of tempMatrix. */
599 int *rowPermutation=(int*)omAlloc(k*sizeof(int));
600 /* This is for storing the permutation of rows
601 resulting from searching for a non-zero
602 pivot element. */
603 for (int i = 0; i < k; i++) rowPermutation[i] = i;
604 int divisor = 1; /* the Bareiss divisor */
605 for (int r = 0; r <= k - 2; r++)
606 {
607 /* look for a non-zero entry in column r: */
608 int i = r;
609 while ((i < k) && (tempMatrix[rowPermutation[i] * k + r] == 0))
610 i++;
611 if (i == k)
612 /* There is no non-zero entry; hence the minor is zero. */
613 return IntMinorValue(0, 0, 0, 0, 0, -1, -1);
614 if (i != r)
615 {
616 /* We swap the rows with indices r and i: */
617 int j = rowPermutation[i];
619 rowPermutation[r] = j;
620 /* Now we know that tempMatrix[rowPermutation[r] * k + r] is not zero.
621 But careful; we have to negate the sign, as there is always an odd
622 number of row transpositions to swap two given rows of a matrix. */
623 sign = -sign;
624 }
625 if (r >= 1) divisor = tempMatrix[rowPermutation[r - 1] * k + r - 1];
626 for (int rr = r + 1; rr < k; rr++)
627 for (int cc = r + 1; cc < k; cc++)
628 {
629 e = rowPermutation[rr] * k + cc;
630 /* Attention: The following may cause an overflow and
631 thus a wrong result: */
633 - tempMatrix[rowPermutation[r] * k + cc]
634 * tempMatrix[rowPermutation[rr] * k + r];
635 /* The following is, by theory, always a division without
636 remainder: */
637 tempMatrix[e] = tempMatrix[e] / divisor;
638 if (characteristic != 0)
640 }
643 }
644 int theValue = sign * tempMatrix[rowPermutation[k - 1] * k + k - 1];
645 if (iSB != 0) theValue = getReduction(theValue, iSB);
646 mv = IntMinorValue(theValue, 0, 0, 0, 0, -1, -1);
647 }
650 return mv;
651}
652
654 const int columnIndex) const
655{
657}
658
660 const int k, const MinorKey& mk,
661 const bool multipleMinors,
663 const int characteristic, const ideal& iSB)
664{
665 assume(k > 0); /* k is the minor's dimension; the minor must be at least
666 1x1 */
667 /* The method works by recursion, and using Lapace's Theorem along
668 the row/column with the most zeros. */
669 if (k == 1)
670 {
671 int e = getEntry(mk.getAbsoluteRowIndex(0), mk.getAbsoluteColumnIndex(0));
672 if (characteristic != 0) e = e % characteristic;
673 if (iSB != 0) e = getReduction(e, iSB);
674 return IntMinorValue(e, 0, 0, 0, 0, -1, -1);
675 /* we set "-1" as, for k == 1, we do not have any cache retrievals */
676 }
677 else
678 {
679 int b = getBestLine(k, mk); /* row or column with
680 most zeros */
681 int result = 0; /* This will contain the
682 value of the minor. */
683 int s = 0; int m = 0; int as = 0; int am = 0; /* counters for additions
684 and multiplications,
685 ..."a*" for
686 accumulated operation
687 counters */
688 IntMinorValue mv(0, 0, 0, 0, 0, 0, 0); /* for storing all
689 intermediate minors */
690 bool hadNonZeroEntry = false;
691 if (b >= 0)
692 {
693 /* This means that the best line is the row with absolute (0-based)
694 index b.
695 Using Laplace, the sign of the contributing minors must be
696 iterating; the initial sign depends on the relative index of b
697 in minorRowKey: */
698 int sign = (mk.getRelativeRowIndex(b) % 2 == 0 ? 1 : -1);
699 for (int c = 0; c < k; c++) /* This iterates over all involved
700 columns. */
701 {
702 int absoluteC = mk.getAbsoluteColumnIndex(c);
703 if (getEntry(b, absoluteC) != 0) /* Only then do we have to consider
704 this sub-determinante. */
705 {
706 hadNonZeroEntry = true;
707 /* Next MinorKey is mk with row b and column absoluteC omitted. */
708 MinorKey subMk = mk.getSubMinorKey(b, absoluteC);
709 if (cch.hasKey(subMk))
710 { /* trying to find the result in the cache */
711 mv = cch.getValue(subMk);
712 mv.incrementRetrievals(); /* once more, we made use of the cached
713 value for key mk */
714 cch.put(subMk, mv); /* We need to do this with "put", as the
715 (altered) number of retrievals may have
716 an impact on the internal ordering among
717 the cached entries. */
718 }
719 else
720 {
722 characteristic, iSB); /* recursive call */
723 /* As this minor was not in the cache, we count the additions
724 and multiplications that we needed to perform in the
725 recursive call: */
726 m += mv.getMultiplications();
727 s += mv.getAdditions();
728 }
729 /* In any case, we count all nested operations in the accumulative
730 counters: */
731 am += mv.getAccumulatedMultiplications();
732 as += mv.getAccumulatedAdditions();
733 /* adding sub-determinante times matrix entry times appropriate
734 sign */
735 result += sign * mv.getResult() * getEntry(b, absoluteC);
737 s++; m++; as++; am++; /* This is for the last addition and
738 multiplication. */
739 }
740 sign = - sign; /* alternating the sign */
741 }
742 }
743 else
744 {
745 b = - b - 1;
746 /* This means that the best line is the column with absolute (0-based)
747 index b.
748 Using Laplace, the sign of the contributing minors must be iterating;
749 the initial sign depends on the relative index of b in
750 minorColumnKey: */
751 int sign = (mk.getRelativeColumnIndex(b) % 2 == 0 ? 1 : -1);
752 for (int r = 0; r < k; r++) /* This iterates over all involved rows. */
753 {
754 int absoluteR = mk.getAbsoluteRowIndex(r);
755 if (getEntry(absoluteR, b) != 0) /* Only then do we have to consider
756 this sub-determinante. */
757 {
758 hadNonZeroEntry = true;
759 /* Next MinorKey is mk with row absoluteR and column b omitted. */
760 MinorKey subMk = mk.getSubMinorKey(absoluteR, b);
761 if (cch.hasKey(subMk))
762 { /* trying to find the result in the cache */
763 mv = cch.getValue(subMk);
764 mv.incrementRetrievals(); /* once more, we made use of the cached
765 value for key mk */
766 cch.put(subMk, mv); /* We need to do this with "put", as the
767 (altered) number of retrievals may have an
768 impact on the internal ordering among the
769 cached entries. */
770 }
771 else
772 {
773 mv = getMinorPrivateLaplace(k - 1, subMk, multipleMinors, cch, characteristic, iSB); /* recursive call */
774 /* As this minor was not in the cache, we count the additions and
775 multiplications that we needed to do in the recursive call: */
776 m += mv.getMultiplications();
777 s += mv.getAdditions();
778 }
779 /* In any case, we count all nested operations in the accumulative
780 counters: */
781 am += mv.getAccumulatedMultiplications();
782 as += mv.getAccumulatedAdditions();
783 /* adding sub-determinante times matrix entry times appropriate
784 sign: */
785 result += sign * mv.getResult() * getEntry(absoluteR, b);
787 s++; m++; as++; am++; /* This is for the last addition and
788 multiplication. */
789 }
790 sign = - sign; /* alternating the sign */
791 }
792 }
793 /* Let's cache the newly computed minor: */
796 _minorSize, k,
798 if (hadNonZeroEntry)
799 {
800 s--; as--; /* first addition was 0 + ..., so we do not count it */
801 }
802 if (s < 0) s = 0; /* may happen when all subminors are zero and no
803 addition needs to be performed */
804 if (as < 0) as = 0; /* may happen when all subminors are zero and no
805 addition needs to be performed */
806 if (iSB != 0) result = getReduction(result, iSB);
808 cch.put(mk, newMV); /* Here's the actual put inside the cache. */
809 return newMV;
810 }
811}
812
817
819 const int columnIndex) const
820{
822}
823
829
831{
832 char h[32];
833 string t = "";
834 string s = "PolyMinorProcessor:";
835 s += "\n matrix: ";
836 sprintf(h, "%d", _rows); s += h;
837 s += " x ";
838 sprintf(h, "%d", _columns); s += h;
839 int myIndexArray[500];
840 s += "\n considered submatrix has row indices: ";
842 for (int k = 0; k < _containerRows; k++)
843 {
844 if (k != 0) s += ", ";
845 sprintf(h, "%d", myIndexArray[k]); s += h;
846 }
847 s += " (first row of matrix has index 0)";
848 s += "\n considered submatrix has column indices: ";
850 for (int k = 0; k < _containerColumns; k++)
851 {
852 if (k != 0) s += ", ";
853 sprintf(h, "%d", myIndexArray[k]); s += h;
854 }
855 s += " (first column of matrix has index 0)";
856 s += "\n size of considered minor(s): ";
857 sprintf(h, "%d", _minorSize); s += h;
858 s += "x";
859 s += h;
860 return s;
861}
862
864{
865 /* free memory of _polyMatrix */
866 int n = _rows * _columns;
867 for (int i = 0; i < n; i++)
870}
871
873 const int numberOfColumns,
874 const poly* polyMatrix)
875{
876 /* free memory of _polyMatrix */
877 int n = _rows * _columns;
878 for (int i = 0; i < n; i++)
881
884 n = _rows * _columns;
885
886 /* allocate memory for new entries in _polyMatrix */
887 _polyMatrix = (poly*)omAlloc(n*sizeof(poly));
888
889 /* copying values from one-dimensional method
890 parameter "polyMatrix" */
891 for (int i = 0; i < n; i++)
893}
894
896 const int* rowIndices,
897 const int* columnIndices,
899 const ideal& iSB)
900{
903 /* call a helper method which recursively computes the minor using the cache
904 c: */
906}
907
909 const int* rowIndices,
910 const int* columnIndices,
911 const char* algorithm,
912 const ideal& iSB)
913{
916 /* call a helper method which computes the minor (without using a cache): */
917 if (strcmp(algorithm, "Laplace") == 0)
919 else if (strcmp(algorithm, "Bareiss") == 0)
921 else assume(false);
922
923 /* The following code is never reached and just there to make the
924 compiler happy: */
925 return PolyMinorValue();
926}
927
929 const ideal& iSB)
930{
931 /* call a helper method which computes the minor (without using a
932 cache): */
933 if (strcmp(algorithm, "Laplace") == 0)
935 else if (strcmp(algorithm, "Bareiss") == 0)
937 else assume(false);
938
939 /* The following code is never reached and just there to make the
940 compiler happy: */
941 return PolyMinorValue();
942}
943
945 PolyMinorValue>& c,
946 const ideal& iSB)
947{
948 /* computation with cache */
949 return getMinorPrivateLaplace(_minorSize, _minor, true, c, iSB);
950}
951
952/* assumes that all entries in polyMatrix are already reduced w.r.t. iSB */
954 const MinorKey& mk,
955 const ideal& iSB)
956{
957 assume(k > 0); /* k is the minor's dimension; the minor must be at least
958 1x1 */
959 /* The method works by recursion, and using Lapace's Theorem along the
960 row/column with the most zeros. */
961 if (k == 1)
962 {
963 PolyMinorValue pmv(getEntry(mk.getAbsoluteRowIndex(0),
964 mk.getAbsoluteColumnIndex(0)),
965 0, 0, 0, 0, -1, -1);
966 /* "-1" is to signal that any statistics about the number of retrievals
967 does not make sense, as we do not use a cache. */
968 return pmv;
969 }
970 else
971 {
972 /* Here, the minor must be 2x2 or larger. */
973 int b = getBestLine(k, mk); /* row or column with most
974 zeros */
975 poly result = NULL; /* This will contain the
976 value of the minor. */
977 int s = 0; int m = 0; int as = 0; int am = 0; /* counters for additions
978 and multiplications,
979 ..."a*" for accumulated
980 operation counters */
981 bool hadNonZeroEntry = false;
982 if (b >= 0)
983 {
984 /* This means that the best line is the row with absolute (0-based)
985 index b.
986 Using Laplace, the sign of the contributing minors must be iterating;
987 the initial sign depends on the relative index of b in minorRowKey: */
988 int sign = (mk.getRelativeRowIndex(b) % 2 == 0 ? 1 : -1);
989 for (int c = 0; c < k; c++) /* This iterates over all involved columns. */
990 {
991 int absoluteC = mk.getAbsoluteColumnIndex(c);
992 if (!isEntryZero(b, absoluteC)) /* Only then do we have to consider
993 this sub-determinante. */
994 {
995 hadNonZeroEntry = true;
996 /* Next MinorKey is mk with row b and column absoluteC omitted. */
997 MinorKey subMk = mk.getSubMinorKey(b, absoluteC);
998 /* recursive call: */
1000 m += mv.getMultiplications();
1001 s += mv.getAdditions();
1002 am += mv.getAccumulatedMultiplications();
1003 as += mv.getAccumulatedAdditions();
1004 poly signPoly = pISet(sign);
1005 poly temp = pp_Mult_qq(mv.getResult(), getEntry(b, absoluteC),
1006 currRing);
1009#ifdef COUNT_AND_PRINT_OPERATIONS
1010 multsPoly++;
1011 addsPoly++;
1012 multsMon += pLength(mv.getResult()) * pLength(getEntry(b, absoluteC));
1013#endif
1014 s++; m++; as++, am++; /* This is for the addition and multiplication
1015 in the previous lines of code. */
1016 }
1017 sign = - sign; /* alternating the sign */
1018 }
1019 }
1020 else
1021 {
1022 b = - b - 1;
1023 /* This means that the best line is the column with absolute (0-based)
1024 index b.
1025 Using Laplace, the sign of the contributing minors must be iterating;
1026 the initial sign depends on the relative index of b in
1027 minorColumnKey: */
1028 int sign = (mk.getRelativeColumnIndex(b) % 2 == 0 ? 1 : -1);
1029 for (int r = 0; r < k; r++) /* This iterates over all involved rows. */
1030 {
1031 int absoluteR = mk.getAbsoluteRowIndex(r);
1032 if (!isEntryZero(absoluteR, b)) /* Only then do we have to consider
1033 this sub-determinante. */
1034 {
1035 hadNonZeroEntry = true;
1036 /* This is mk with row absoluteR and column b omitted. */
1037 MinorKey subMk = mk.getSubMinorKey(absoluteR, b);
1038 /* recursive call: */
1040 m += mv.getMultiplications();
1041 s += mv.getAdditions();
1042 am += mv.getAccumulatedMultiplications();
1043 as += mv.getAccumulatedAdditions();
1044 poly signPoly = pISet(sign);
1045 poly temp = pp_Mult_qq(mv.getResult(), getEntry(absoluteR, b),
1046 currRing);
1049#ifdef COUNT_AND_PRINT_OPERATIONS
1050 multsPoly++;
1051 addsPoly++;
1052 multsMon += pLength(mv.getResult()) * pLength(getEntry(absoluteR, b));
1053#endif
1054 s++; m++; as++, am++; /* This is for the addition and multiplication
1055 in the previous lines of code. */
1056 }
1057 sign = - sign; /* alternating the sign */
1058 }
1059 }
1060 if (hadNonZeroEntry)
1061 {
1062 s--; as--; /* first addition was 0 + ..., so we do not count it */
1063#ifdef COUNT_AND_PRINT_OPERATIONS
1064 addsPoly--;
1065#endif
1066 }
1067 if (s < 0) s = 0; /* may happen when all subminors are zero and no
1068 addition needs to be performed */
1069 if (as < 0) as = 0; /* may happen when all subminors are zero and no
1070 addition needs to be performed */
1071 if (iSB != NULL)
1072 {
1073 poly tmpresult = kNF(iSB, currRing->qideal, result);
1074 pDelete(&result);
1076 }
1077 PolyMinorValue newMV(result, m, s, am, as, -1, -1);
1078 /* "-1" is to signal that any statistics about the number of retrievals
1079 does not make sense, as we do not use a cache. */
1080 pDelete(&result);
1081 return newMV;
1082 }
1083}
1084
1085/* assumes that all entries in polyMatrix are already reduced w.r.t. iSB */
1087 const int k,
1088 const MinorKey& mk,
1089 const bool multipleMinors,
1091 const ideal& iSB)
1092{
1093 assume(k > 0); /* k is the minor's dimension; the minor must be at least
1094 1x1 */
1095 /* The method works by recursion, and using Lapace's Theorem along
1096 the row/column with the most zeros. */
1097 if (k == 1)
1098 {
1099 PolyMinorValue pmv(getEntry(mk.getAbsoluteRowIndex(0),
1100 mk.getAbsoluteColumnIndex(0)),
1101 0, 0, 0, 0, -1, -1);
1102 /* we set "-1" as, for k == 1, we do not have any cache retrievals */
1103 return pmv;
1104 }
1105 else
1106 {
1107 int b = getBestLine(k, mk); /* row or column with most
1108 zeros */
1109 poly result = NULL; /* This will contain the
1110 value of the minor. */
1111 int s = 0; int m = 0; int as = 0; int am = 0; /* counters for additions
1112 and multiplications,
1113 ..."a*" for accumulated
1114 operation counters */
1115 bool hadNonZeroEntry = false;
1116 if (b >= 0)
1117 {
1118 /* This means that the best line is the row with absolute (0-based)
1119 index b.
1120 Using Laplace, the sign of the contributing minors must be iterating;
1121 the initial sign depends on the relative index of b in
1122 minorRowKey: */
1123 int sign = (mk.getRelativeRowIndex(b) % 2 == 0 ? 1 : -1);
1124 for (int c = 0; c < k; c++) /* This iterates over all involved columns. */
1125 {
1126 int absoluteC = mk.getAbsoluteColumnIndex(c);
1127 if (!isEntryZero(b, absoluteC)) /* Only then do we have to consider
1128 this sub-determinante. */
1129 {
1130 hadNonZeroEntry = true;
1131 PolyMinorValue mv; /* for storing all intermediate minors */
1132 /* Next MinorKey is mk with row b and column absoluteC omitted. */
1133 MinorKey subMk = mk.getSubMinorKey(b, absoluteC);
1134 if (cch.hasKey(subMk))
1135 { /* trying to find the result in the cache */
1136 mv = cch.getValue(subMk);
1137 mv.incrementRetrievals(); /* once more, we made use of the cached
1138 value for key mk */
1139 cch.put(subMk, mv); /* We need to do this with "put", as the
1140 (altered) number of retrievals may have an
1141 impact on the internal ordering among cache
1142 entries. */
1143 }
1144 else
1145 {
1146 /* recursive call: */
1148 iSB);
1149 /* As this minor was not in the cache, we count the additions and
1150 multiplications that we needed to do in the recursive call: */
1151 m += mv.getMultiplications();
1152 s += mv.getAdditions();
1153 }
1154 /* In any case, we count all nested operations in the accumulative
1155 counters: */
1156 am += mv.getAccumulatedMultiplications();
1157 as += mv.getAccumulatedAdditions();
1158 poly signPoly = pISet(sign);
1159 poly temp = pp_Mult_qq(mv.getResult(), getEntry(b, absoluteC),
1160 currRing);
1163#ifdef COUNT_AND_PRINT_OPERATIONS
1164 multsPoly++;
1165 addsPoly++;
1166 multsMon += pLength(mv.getResult()) * pLength(getEntry(b, absoluteC));
1167#endif
1168 s++; m++; as++; am++; /* This is for the addition and multiplication
1169 in the previous lines of code. */
1170 }
1171 sign = - sign; /* alternating the sign */
1172 }
1173 }
1174 else
1175 {
1176 b = - b - 1;
1177 /* This means that the best line is the column with absolute (0-based)
1178 index b.
1179 Using Laplace, the sign of the contributing minors must be iterating;
1180 the initial sign depends on the relative index of b in
1181 minorColumnKey: */
1182 int sign = (mk.getRelativeColumnIndex(b) % 2 == 0 ? 1 : -1);
1183 for (int r = 0; r < k; r++) /* This iterates over all involved rows. */
1184 {
1185 int absoluteR = mk.getAbsoluteRowIndex(r);
1186 if (!isEntryZero(absoluteR, b)) /* Only then do we have to consider
1187 this sub-determinante. */
1188 {
1189 hadNonZeroEntry = true;
1190 PolyMinorValue mv; /* for storing all intermediate minors */
1191 /* Next MinorKey is mk with row absoluteR and column b omitted. */
1192 MinorKey subMk = mk.getSubMinorKey(absoluteR, b);
1193 if (cch.hasKey(subMk))
1194 { /* trying to find the result in the cache */
1195 mv = cch.getValue(subMk);
1196 mv.incrementRetrievals(); /* once more, we made use of the cached
1197 value for key mk */
1198 cch.put(subMk, mv); /* We need to do this with "put", as the
1199 (altered) number of retrievals may have an
1200 impact on the internal ordering among the
1201 cached entries. */
1202 }
1203 else
1204 {
1206 iSB); /* recursive call */
1207 /* As this minor was not in the cache, we count the additions and
1208 multiplications that we needed to do in the recursive call: */
1209 m += mv.getMultiplications();
1210 s += mv.getAdditions();
1211 }
1212 /* In any case, we count all nested operations in the accumulative
1213 counters: */
1214 am += mv.getAccumulatedMultiplications();
1215 as += mv.getAccumulatedAdditions();
1216 poly signPoly = pISet(sign);
1217 poly temp = pp_Mult_qq(mv.getResult(), getEntry(absoluteR, b),
1218 currRing);
1221#ifdef COUNT_AND_PRINT_OPERATIONS
1222 multsPoly++;
1223 addsPoly++;
1224 multsMon += pLength(mv.getResult()) * pLength(getEntry(absoluteR, b));
1225#endif
1226 s++; m++; as++; am++; /* This is for the addition and multiplication
1227 in the previous lines of code. */
1228 }
1229 sign = - sign; /* alternating the sign */
1230 }
1231 }
1232 /* Let's cache the newly computed minor: */
1235 _minorSize,
1236 k,
1238 if (hadNonZeroEntry)
1239 {
1240 s--; as--; /* first addition was 0 + ..., so we do not count it */
1241#ifdef COUNT_AND_PRINT_OPERATIONS
1242 addsPoly--;
1243#endif
1244 }
1245 if (s < 0) s = 0; /* may happen when all subminors are zero and no
1246 addition needs to be performed */
1247 if (as < 0) as = 0; /* may happen when all subminors are zero and no
1248 addition needs to be performed */
1249 if (iSB != NULL)
1250 {
1251 poly tmpresult = kNF(iSB, currRing->qideal, result);
1254 }
1257 cch.put(mk, newMV); /* Here's the actual put inside the cache. */
1258 return newMV;
1259 }
1260}
1261
1262/* This can only be used in the case of coefficients coming from a field
1263 or at least an integral domain. */
1264static void addOperationBucket(poly f1, poly f2, kBucket_pt bucket)
1265{
1266 /* fills all terms of f1 * f2 into the bucket */
1267 poly a = f1; poly b = f2;
1268 int aLen = pLength(a); int bLen = pLength(b);
1269 if (aLen > bLen)
1270 {
1271 b = f1; a = f2; bLen = aLen;
1272 }
1273 pNormalize(b);
1274
1275 while (a != NULL)
1276 {
1277 /* The next line actually uses only LT(a): */
1278 kBucket_Plus_mm_Mult_pp(bucket, a, b, bLen);
1279 a = pNext(a);
1280 }
1281}
1282
1283/* computes the polynomial (p1 * p2 - p3 * p4) and puts result into p1;
1284 the method destroys the old value of p1;
1285 p2, p3, and p4 may be pNormalize-d but must, apart from that,
1286 not be changed;
1287 This can only be used in the case of coefficients coming from a field
1288 or at least an integral domain. */
1289static void elimOperationBucketNoDiv(poly &p1, poly p2, poly p3, poly p4)
1290{
1291#ifdef COUNT_AND_PRINT_OPERATIONS
1292 if ((pLength(p1) != 0) && (pLength(p2) != 0))
1293 {
1294 multsPoly++;
1295 multsMon += pLength(p1) * pLength(p2);
1296 }
1297 if ((pLength(p3) != 0) && (pLength(p4) != 0))
1298 {
1299 multsPoly++;
1300 multsMon += pLength(p3) * pLength(p4);
1301 }
1302 if ((pLength(p1) != 0) && (pLength(p2) != 0) &&
1303 (pLength(p3) != 0) && (pLength(p4) != 0))
1304 addsPoly++;
1305#endif
1308 poly p3Neg = pNeg(pCopy(p3));
1310 pDelete(&p3Neg);
1311 pDelete(&p1);
1312 p1 = kBucketClear(myBucket);
1314}
1315
1316/* computes the polynomial (p1 * p2 - p3 * p4) / p5 and puts result into p1;
1317 the method destroys the old value of p1;
1318 p2, p3, p4, and p5 may be pNormalize-d but must, apart from that,
1319 not be changed;
1320 c5 is assumed to be the leading coefficient of p5;
1321 p5Len is assumed to be the length of p5;
1322 This can only be used in the case of coefficients coming from a field
1323 or at least an integral domain. */
1324void elimOperationBucket(poly &p1, poly &p2, poly &p3, poly &p4, poly &p5,
1325 number &c5, int p5Len)
1326{
1327#ifdef COUNT_AND_PRINT_OPERATIONS
1328 if ((pLength(p1) != 0) && (pLength(p2) != 0))
1329 {
1330 multsPoly++;
1331 multsMon += pLength(p1) * pLength(p2);
1332 }
1333 if ((pLength(p3) != 0) && (pLength(p4) != 0))
1334 {
1335 multsPoly++;
1336 multsMon += pLength(p3) * pLength(p4);
1337 }
1338 if ((pLength(p1) != 0) && (pLength(p2) != 0) &&
1339 (pLength(p3) != 0) && (pLength(p4) != 0))
1340 addsPoly++;
1341#endif
1344 poly p3Neg = pNeg(pCopy(p3));
1346 pDelete(&p3Neg);
1347
1348 /* Now, myBucket contains all terms of p1 * p2 - p3 * p4.
1349 Now we need to perform the polynomial division myBucket / p5
1350 which is known to work without remainder: */
1351 pDelete(&p1); poly helperPoly = NULL;
1352
1354 while (bucketLm != NULL)
1355 {
1356 /* divide bucketLm by the leading term of p5 and put result into bucketLm;
1357 we start with the coefficients;
1358 note that bucketLm will always represent a term */
1359 number coeff = nDiv(pGetCoeff(bucketLm), c5);
1360 nNormalize(coeff);
1361 pSetCoeff(bucketLm, coeff);
1362 /* subtract exponent vector of p5 from that of quotient; modifies
1363 quotient */
1365#ifdef COUNT_AND_PRINT_OPERATIONS
1366 divsMon++;
1368 multsMon += p5Len;
1369 savedMultsMFD++;
1370 multsPoly++;
1372 addsPoly++;
1374#endif
1376 /* The following lines make bucketLm the new leading term of p1,
1377 i.e., put bucketLm in front of everything which is already in p1.
1378 Thus, after the while loop, we need to revert p1. */
1380 helperPoly->next = p1;
1381 p1 = helperPoly;
1382
1384 }
1385 p1 = pReverse(p1);
1387}
1388
1389/* assumes that all entries in polyMatrix are already reduced w.r.t. iSB
1390 This can only be used in the case of coefficients coming from a field!!! */
1392 const MinorKey& mk,
1393 const ideal& iSB)
1394{
1395 assume(k > 0); /* k is the minor's dimension; the minor must be at least
1396 1x1 */
1397 int *theRows=(int*)omAlloc(k*sizeof(int));
1398 mk.getAbsoluteRowIndices(theRows);
1399 int *theColumns=(int*)omAlloc(k*sizeof(int));
1400 mk.getAbsoluteColumnIndices(theColumns);
1401 if (k == 1)
1402 {
1404 0, 0, 0, 0, -1, -1);
1406 omFree(theRows);
1407 return tmp;
1408 }
1409 else /* k > 0 */
1410 {
1411 /* the matrix to perform Bareiss with */
1412 poly* tempMatrix = (poly*)omAlloc(k * k * sizeof(poly));
1413 /* copy correct set of entries from _polyMatrix to tempMatrix */
1414 int i = 0;
1415 for (int r = 0; r < k; r++)
1416 for (int c = 0; c < k; c++)
1418
1419 /* Bareiss algorithm operating on tempMatrix which is at least 2x2 */
1420 int sign = 1; /* This will store the correct sign resulting from
1421 permuting the rows of tempMatrix. */
1422 int *rowPermutation=(int*)omAlloc(k*sizeof(int));
1423 /* This is for storing the permutation of rows
1424 resulting from searching for a non-zero pivot
1425 element. */
1426 for (int i = 0; i < k; i++) rowPermutation[i] = i;
1427 poly divisor = NULL;
1428 int divisorLength = 0;
1430 for (int r = 0; r <= k - 2; r++)
1431 {
1432 /* look for a non-zero entry in column r, rows = r .. (k - 1)
1433 s.t. the polynomial has least complexity: */
1434 int minComplexity = -1; int complexity = 0; int bestRow = -1;
1435 poly pp = NULL;
1436 for (int i = r; i < k; i++)
1437 {
1438 pp = tempMatrix[rowPermutation[i] * k + r];
1439 if (pp != NULL)
1440 {
1441 if (minComplexity == -1)
1442 {
1444 bestRow = i;
1445 }
1446 else
1447 {
1448 complexity = 0;
1449 while ((pp != NULL) && (complexity < minComplexity))
1450 {
1451 complexity += nSize(pGetCoeff(pp)); pp = pNext(pp);
1452 }
1453 if (complexity < minComplexity)
1454 {
1455 minComplexity = complexity;
1456 bestRow = i;
1457 }
1458 }
1459 if (minComplexity <= 1) break; /* terminate the search */
1460 }
1461 }
1462 if (bestRow == -1)
1463 {
1464 /* There is no non-zero entry; hence the minor is zero. */
1465 for (int i = 0; i < k * k; i++) pDelete(&tempMatrix[i]);
1466 return PolyMinorValue(NULL, 0, 0, 0, 0, -1, -1);
1467 }
1469 if (bestRow != r)
1470 {
1471 /* We swap the rows with indices r and i: */
1472 int j = rowPermutation[bestRow];
1474 rowPermutation[r] = j;
1475 /* Now we know that tempMatrix[rowPermutation[r] * k + r] is not zero.
1476 But careful; we have to negate the sign, as there is always an odd
1477 number of row transpositions to swap two given rows of a matrix. */
1478 sign = -sign;
1479 }
1480#if (defined COUNT_AND_PRINT_OPERATIONS) && (COUNT_AND_PRINT_OPERATIONS > 2)
1481 poly w = NULL; int wl = 0;
1482 printf("matrix after %d steps:\n", r);
1483 for (int u = 0; u < k; u++)
1484 {
1485 for (int v = 0; v < k; v++)
1486 {
1487 if ((v < r) && (u > v))
1488 wl = 0;
1489 else
1490 {
1491 w = tempMatrix[rowPermutation[u] * k + v];
1492 wl = pLength(w);
1493 }
1494 printf("%5d ", wl);
1495 }
1496 printf("\n");
1497 }
1498 printCounters ("", false);
1499#endif
1500 if (r != 0)
1501 {
1502 divisor = tempMatrix[rowPermutation[r - 1] * k + r - 1];
1503 pNormalize(divisor);
1504 divisorLength = pLength(divisor);
1505 divisorLC = pGetCoeff(divisor);
1506 }
1507 for (int rr = r + 1; rr < k; rr++)
1508 for (int cc = r + 1; cc < k; cc++)
1509 {
1510 if (r == 0)
1512 tempMatrix[rowPermutation[r] * k + r],
1513 tempMatrix[rowPermutation[r] * k + cc],
1514 tempMatrix[rowPermutation[rr] * k + r]);
1515 else
1517 tempMatrix[rowPermutation[r] * k + r],
1518 tempMatrix[rowPermutation[r] * k + cc],
1519 tempMatrix[rowPermutation[rr] * k + r],
1520 divisor, divisorLC, divisorLength);
1521 }
1522 }
1523#if (defined COUNT_AND_PRINT_OPERATIONS) && (COUNT_AND_PRINT_OPERATIONS > 2)
1524 poly w = NULL; int wl = 0;
1525 printf("matrix after %d steps:\n", k - 1);
1526 for (int u = 0; u < k; u++)
1527 {
1528 for (int v = 0; v < k; v++)
1529 {
1530 if ((v < k - 1) && (u > v))
1531 wl = 0;
1532 else
1533 {
1534 w = tempMatrix[rowPermutation[u] * k + v];
1535 wl = pLength(w);
1536 }
1537 printf("%5d ", wl);
1538 }
1539 printf("\n");
1540 }
1541#endif
1542 poly result = tempMatrix[rowPermutation[k - 1] * k + k - 1];
1543 tempMatrix[rowPermutation[k - 1] * k + k - 1]=NULL; /*value now in result*/
1544 if (sign == -1) result = pNeg(result);
1545 if (iSB != NULL)
1546 {
1547 poly tmpresult = kNF(iSB, currRing->qideal, result);
1548 pDelete(&result);
1550 }
1551 PolyMinorValue mv(result, 0, 0, 0, 0, -1, -1);
1552 for (int i = 0; i < k * k; i++) pDelete(&tempMatrix[i]);
1553 omFreeSize(tempMatrix, k * k * sizeof(poly));
1556 omFree(theRows);
1557 return mv;
1558 }
1559}
1560
static void addOperationBucket(poly f1, poly f2, kBucket_pt bucket)
int getReduction(const int i, const ideal &iSB)
static void elimOperationBucketNoDiv(poly &p1, poly p2, poly p3, poly p4)
void elimOperationBucket(poly &p1, poly &p2, poly &p3, poly &p4, poly &p5, number &c5, int p5Len)
void printCounters(char *prefix, bool resetToZero)
BOOLEAN dimension(leftv res, leftv args)
Definition bbcone.cc:757
CanonicalForm FACTORY_PUBLIC pp(const CanonicalForm &)
CanonicalForm pp ( const CanonicalForm & f )
Definition cf_gcd.cc:676
int m
Definition cfEzgcd.cc:128
int i
Definition cfEzgcd.cc:132
int k
Definition cfEzgcd.cc:99
g
Definition cfModGcd.cc:4091
CanonicalForm b
Definition cfModGcd.cc:4104
FILE * f
Definition checklibs.c:9
Class Cache is a template-implementation of a cache with arbitrary classes for representing keys and ...
Definition Cache.h:69
std::string toString() const
A method for providing a printable version of the represented MinorProcessor.
~IntMinorProcessor()
A destructor for deleting an instance.
IntMinorProcessor()
A constructor for creating an instance.
int * _intMatrix
private store for integer matrix entries
IntMinorValue getMinorPrivateBareiss(const int k, const MinorKey &mk, const int characteristic, const ideal &iSB)
A method for computing the value of a minor using Bareiss's algorithm.
bool isEntryZero(const int absoluteRowIndex, const int absoluteColumnIndex) const
A method for testing whether a matrix entry is zero.
void defineMatrix(const int numberOfRows, const int numberOfColumns, const int *matrix)
A method for defining a matrix with integer entries.
IntMinorValue getMinor(const int dimension, const int *rowIndices, const int *columnIndices, const int characteristic, const ideal &iSB, const char *algorithm)
A method for computing the value of a minor without using a cache.
IntMinorValue getNextMinor(const int characteristic, const ideal &iSB, const char *algorithm)
A method for obtaining the next minor when iterating through all minors of a given size within a pre-...
IntMinorValue getMinorPrivateLaplace(const int k, const MinorKey &mk, const bool multipleMinors, Cache< MinorKey, IntMinorValue > &c, int characteristic, const ideal &iSB)
A method for computing the value of a minor, using a cache.
int getEntry(const int rowIndex, const int columnIndex) const
A method for retrieving the matrix entry.
Class IntMinorValue is derived from MinorValue and can be used for representing values in a cache for...
Definition Minor.h:718
Class MinorKey can be used for representing keys in a cache for sub-determinantes; see class Cache.
Definition Minor.h:40
void getAbsoluteColumnIndices(int *const target) const
A method for retrieving the 0-based indices of all columns encoded in this MinorKey.
Definition Minor.cc:202
void selectFirstRows(const int k, const MinorKey &mk)
This method redefines the set of rows represented by this MinorKey.
Definition Minor.cc:457
bool selectNextColumns(const int k, const MinorKey &mk)
This method redefines the set of columns represented by this MinorKey.
Definition Minor.cc:669
void reset()
A method for deleting all entries of _rowKey and _columnKey.
Definition Minor.cc:13
void selectFirstColumns(const int k, const MinorKey &mk)
This method redefines the set of columns represented by this MinorKey.
Definition Minor.cc:498
void set(const int lengthOfRowArray, const unsigned int *rowKey, const int lengthOfColumnArray, const unsigned int *columnKey)
A setter method for class MinorKey.
Definition Minor.cc:62
void getAbsoluteRowIndices(int *const target) const
A method for retrieving the 0-based indices of all rows encoded in this MinorKey.
Definition Minor.cc:181
bool selectNextRows(const int k, const MinorKey &mk)
This method redefines the set of rows represented by this MinorKey.
Definition Minor.cc:538
int compare(const MinorKey &mk) const
A comparator for two instances of MinorKey.
Definition Minor.cc:412
virtual bool isEntryZero(const int absoluteRowIndex, const int absoluteColumnIndex) const
A method for testing whether a matrix entry is zero.
static int IOverJ(const int i, const int j)
A static method for computing the binomial coefficient i over j.
MinorKey _minor
private store for the rows and columns of the minor of interest; Usually, this minor will encode subs...
void getCurrentColumnIndices(int *const target) const
A method for obtaining the current set of columns corresponding to the current minor when iterating t...
static int NumberOfRetrievals(const int rows, const int columns, const int containerMinorSize, const int minorSize, const bool multipleMinors)
A static method for computing the maximum number of retrievals of a minor.
static int Faculty(const int i)
A static method for computing the factorial of i.
void setMinorSize(const int minorSize)
Sets the size of the minor(s) of interest.
int _containerRows
private store for the number of rows in the container minor; This is set by MinorProcessor::defineSub...
virtual std::string toString() const
A method for providing a printable version of the represented MinorProcessor.
int getBestLine(const int k, const MinorKey &mk) const
A method for identifying the row or column with the most zeros.
bool setNextKeys(const int k)
A method for iterating through all possible subsets of k rows and k columns inside a pre-defined subm...
void print() const
A method for printing a string representation of the given MinorProcessor to std::cout.
int _columns
private store for the number of columns in the underlying matrix
int _minorSize
private store for the dimension of the minor(s) of interest
void defineSubMatrix(const int numberOfRows, const int *rowIndices, const int numberOfColumns, const int *columnIndices)
A method for defining a sub-matrix within a pre-defined matrix.
MinorKey _container
private store for the rows and columns of the container minor within the underlying matrix; _containe...
bool hasNextMinor()
A method for checking whether there is a next choice of rows and columns when iterating through all m...
virtual ~MinorProcessor()
A destructor for deleting an instance.
void getCurrentRowIndices(int *const target) const
A method for obtaining the current set of rows corresponding to the current minor when iterating thro...
MinorProcessor()
The default constructor.
int _rows
private store for the number of rows in the underlying matrix
int _containerColumns
private store for the number of columns in the container minor; This is set by MinorProcessor::define...
~PolyMinorProcessor()
A destructor for deleting an instance.
std::string toString() const
A method for providing a printable version of the represented MinorProcessor.
PolyMinorProcessor()
A constructor for creating an instance.
PolyMinorValue getNextMinor(const char *algorithm, const ideal &iSB)
A method for obtaining the next minor when iterating through all minors of a given size within a pre-...
bool isEntryZero(const int absoluteRowIndex, const int absoluteColumnIndex) const
A method for testing whether a matrix entry is zero.
poly getEntry(const int rowIndex, const int columnIndex) const
A method for retrieving the matrix entry.
void defineMatrix(const int numberOfRows, const int numberOfColumns, const poly *polyMatrix)
A method for defining a matrix with polynomial entries.
PolyMinorValue getMinor(const int dimension, const int *rowIndices, const int *columnIndices, const char *algorithm, const ideal &iSB)
A method for computing the value of a minor, without using a cache.
PolyMinorValue getMinorPrivateLaplace(const int k, const MinorKey &mk, const bool multipleMinors, Cache< MinorKey, PolyMinorValue > &c, const ideal &iSB)
A method for computing the value of a minor, using a cache.
PolyMinorValue getMinorPrivateBareiss(const int k, const MinorKey &mk, const ideal &iSB)
A method for computing the value of a minor, without using a cache.
poly * _polyMatrix
private store for polynomial matrix entries
Class PolyMinorValue is derived from MinorValue and can be used for representing values in a cache fo...
Definition Minor.h:800
static FORCE_INLINE long n_Int(number &n, const coeffs r)
conversion of n to an int; 0 if not possible in Z/pZ: the representing int lying in (-p/2 ....
Definition coeffs.h:544
return result
const CanonicalForm int s
Definition facAbsFact.cc:51
const CanonicalForm & w
Definition facAbsFact.cc:51
const Variable & v
< [in] a sqrfree bivariate poly
Definition facBivar.h:39
int j
Definition facHensel.cc:110
#define VAR
Definition globaldefs.h:5
int binom(int n, int r)
STATIC_VAR int offset
Definition janet.cc:29
STATIC_VAR Poly * h
Definition janet.cc:971
void kBucketClear(kBucket_pt bucket, poly *p, int *length)
Definition kbuckets.cc:521
void kBucket_Minus_m_Mult_p(kBucket_pt bucket, poly m, poly p, int *l, poly spNoether)
Bpoly == Bpoly - m*p; where m is a monom Does not destroy p and m assume (*l <= 0 || pLength(p) == *l...
Definition kbuckets.cc:722
void kBucketDestroy(kBucket_pt *bucket_pt)
Definition kbuckets.cc:216
kBucket_pt kBucketCreate(const ring bucket_ring)
Creation/Destruction of buckets.
Definition kbuckets.cc:209
void kBucket_Plus_mm_Mult_pp(kBucket_pt bucket, poly m, poly p, int l)
Bpoly == Bpoly + m*p; where m is a monom Does not destroy p and m assume (l <= 0 || pLength(p) == l)
Definition kbuckets.cc:815
const poly kBucketGetLm(kBucket_pt bucket)
Definition kbuckets.cc:506
poly kNF(ideal F, ideal Q, poly p, int syzComp, int lazyReduce)
Definition kstd1.cc:3235
#define assume(x)
Definition mod2.h:389
#define pNext(p)
Definition monomials.h:36
static number & pGetCoeff(poly p)
return an alias to the leading coefficient of p assumes that p != NULL NOTE: not copy
Definition monomials.h:44
#define nDiv(a, b)
Definition numbers.h:32
#define nSize(n)
Definition numbers.h:39
#define nNormalize(n)
Definition numbers.h:30
#define omfree(addr)
#define omFreeSize(addr, size)
#define omAlloc(size)
#define omFree(addr)
#define omAlloc0(size)
#define NULL
Definition omList.c:12
static int pLength(poly a)
Definition p_polys.h:190
static poly p_Add_q(poly p, poly q, const ring r)
Definition p_polys.h:936
static poly p_Mult_q(poly p, poly q, const ring r)
Definition p_polys.h:1114
static void p_ExpVectorSub(poly p1, poly p2, const ring r)
Definition p_polys.h:1440
static poly pReverse(poly p)
Definition p_polys.h:335
static void p_Delete(poly *p, const ring r)
Definition p_polys.h:901
static poly pp_Mult_qq(poly p, poly q, const ring r)
Definition p_polys.h:1151
VAR ring currRing
Widely used global variable which specifies the current polynomial ring for Singular interpreter and ...
Definition polys.cc:13
Compatibility layer for legacy polynomial operations (over currRing)
#define pDelete(p_ptr)
Definition polys.h:186
#define pNeg(p)
Definition polys.h:198
#define pSetCoeff(p, n)
deletes old coeff before setting the new one
Definition polys.h:31
#define pNormalize(p)
Definition polys.h:317
#define pSize(p)
Definition polys.h:318
#define pCopy(p)
return a copy of the poly
Definition polys.h:185
#define pISet(i)
Definition polys.h:312
void PrintS(const char *s)
Definition reporter.cc:284
static int sign(int x)
Definition ring.cc:3442