18#ifdef COUNT_AND_PRINT_OPERATIONS
20VAR long multsPoly = 0;
21VAR long addsPolyForDiv = 0;
23VAR long multsPolyForDiv = 0;
26VAR long multsMonForDiv = 0;
28VAR long savedMultsMFD = 0;
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",
40 addsPoly, addsPolyForDiv, multsPoly, multsPolyForDiv,
41 multsMon, multsMonForDiv, savedMultsMFD, divsMon);
44 multsMon = 0; addsPoly = 0; multsPoly = 0; divsMon = 0;
45 savedMultsMFD = 0; multsMonForDiv = 0; addsPolyForDiv = 0;
66 int numberOfZeros = 0;
67 int bestIndex = 100000;
68 int maxNumberOfZeros = -1;
70 for (
int r = 0; r <
k; r++)
75 for (
int c = 0; c <
k; c++)
78 if (
isEntryZero(absoluteR, absoluteC)) numberOfZeros++;
80 if (numberOfZeros > maxNumberOfZeros)
83 bestIndex = absoluteR;
84 maxNumberOfZeros = numberOfZeros;
87 for (
int c = 0; c <
k; c++)
91 for (
int r = 0; r <
k; r++)
94 if (
isEntryZero(absoluteR, absoluteC)) numberOfZeros++;
96 if (numberOfZeros > maxNumberOfZeros)
101 bestIndex = - absoluteC - 1;
102 maxNumberOfZeros = numberOfZeros;
130 const int* rowIndices,
131 const int numberOfColumns,
132 const int* columnIndices)
142 int highestRowIndex = rowIndices[numberOfRows - 1];
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++)
148 int blockIndex = rowIndices[
i] / 32;
149 int offset = rowIndices[
i] % 32;
150 rowBlocks[blockIndex] += (1 <<
offset);
154 int highestColumnIndex = columnIndices[numberOfColumns - 1];
155 int columnBlockCount = (highestColumnIndex / 32) + 1;
156 unsigned *columnBlocks=(
unsigned*)
omAlloc0(columnBlockCount*
sizeof(
unsigned));
157 for (
int i = 0;
i < numberOfColumns;
i++)
159 int blockIndex = columnIndices[
i] / 32;
160 int offset = columnIndices[
i] % 32;
161 columnBlocks[blockIndex] += (1 <<
offset);
164 _container.
set(rowBlockCount, rowBlocks, columnBlockCount, columnBlocks);
222 if (
j == 0 ||
i ==
j)
return 1;
246 const int containerMinorSize,
248 const bool multipleMinors)
262 result =
IOverJ(rows - minorSize, containerMinorSize - minorSize)
263 *
IOverJ(columns - minorSize, containerMinorSize - minorSize)
264 *
Faculty(containerMinorSize - minorSize);
281 _containerColumns(0),
299 string s =
"IntMinorProcessor:";
304 for (
int r = 0; r <
_rows; r++)
310 for (
int k = 0;
k < int(4 - strlen(
h));
k++)
s +=
" ";
314 int myIndexArray[500];
315 s +=
"\n considered submatrix has row indices: ";
319 if (
k != 0)
s +=
", ";
320 sprintf(
h,
"%d", myIndexArray[
k]);
s +=
h;
322 s +=
" (first row of matrix has index 0)";
323 s +=
"\n considered submatrix has column indices: ";
327 if (
k != 0)
s +=
", ";
328 sprintf(
h,
"%d", myIndexArray[
k]);
s +=
h;
330 s +=
" (first column of matrix has index 0)";
331 s +=
"\n size of considered minor(s): ";
345 const int absoluteColumnIndex)
const
347 return getEntry(absoluteRowIndex, absoluteColumnIndex) == 0;
351 const int numberOfColumns,
357 _rows = numberOfRows;
366 for (
int i = 0;
i < n;
i++)
371 const int* rowIndices,
372 const int* columnIndices,
374 const int characteristic,
382 characteristic, iSB);
386 const int* rowIndices,
387 const int* columnIndices,
388 const int characteristic,
390 const char* algorithm)
396 if (strcmp(algorithm,
"Laplace") == 0)
399 else if (strcmp(algorithm,
"Bareiss") == 0)
411 const char* algorithm)
414 if (strcmp(algorithm,
"Laplace") == 0)
416 else if (strcmp(algorithm,
"Bareiss") == 0)
427 const int characteristic,
439 if (
i == 0)
return 0;
452 const int characteristic,
462 if (characteristic != 0) e = e % characteristic;
477 int s = 0;
int m = 0;
int as = 0;
int am = 0;
481 bool hadNonZeroEntry =
false;
489 for (
int c = 0; c <
k; c++)
495 hadNonZeroEntry =
true;
499 characteristic, iSB);
508 if (characteristic != 0)
result =
result % characteristic;
509 s++;
m++; as++, am++;
524 for (
int r = 0; r <
k; r++)
530 hadNonZeroEntry =
true;
541 if (characteristic != 0)
result =
result % characteristic;
542 s++;
m++; as++, am++;
569 const int characteristic,
574 int *theRows=(
int*)
omAlloc(
k*
sizeof(
int));
576 int *theColumns=(
int*)
omAlloc(
k*
sizeof(
int));
579 int e =
getEntry(theRows[0], theColumns[0]);
580 if (characteristic != 0) e = e % characteristic;
586 long *tempMatrix=(
long*)
omAlloc(
k *
k *
sizeof(
long));
589 for (
int r = 0; r <
k; r++)
590 for (
int c = 0; c <
k; c++)
592 e =
getEntry(theRows[r], theColumns[c]);
593 if (characteristic != 0) e = e % characteristic;
599 int *rowPermutation=(
int*)
omAlloc(
k*
sizeof(
int));
603 for (
int i = 0;
i <
k;
i++) rowPermutation[
i] =
i;
605 for (
int r = 0; r <=
k - 2; r++)
609 while ((
i <
k) && (tempMatrix[rowPermutation[
i] *
k + r] == 0))
617 int j = rowPermutation[
i];
618 rowPermutation[
i] = rowPermutation[r];
619 rowPermutation[r] =
j;
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++)
629 e = rowPermutation[rr] *
k + cc;
632 tempMatrix[e] = tempMatrix[e] * tempMatrix[rowPermutation[r] *
k + r]
633 - tempMatrix[rowPermutation[r] *
k + cc]
634 * tempMatrix[rowPermutation[rr] *
k + r];
637 tempMatrix[e] = tempMatrix[e] / divisor;
638 if (characteristic != 0)
639 tempMatrix[e] = tempMatrix[e] % characteristic;
644 int theValue =
sign * tempMatrix[rowPermutation[
k - 1] *
k +
k - 1];
654 const int columnIndex)
const
661 const bool multipleMinors,
663 const int characteristic,
const ideal& iSB)
672 if (characteristic != 0) e = e % characteristic;
683 int s = 0;
int m = 0;
int as = 0;
int am = 0;
690 bool hadNonZeroEntry =
false;
699 for (
int c = 0; c <
k; c++)
706 hadNonZeroEntry =
true;
722 characteristic, iSB);
736 if (characteristic != 0)
result =
result % characteristic;
737 s++;
m++; as++; am++;
752 for (
int r = 0; r <
k; r++)
758 hadNonZeroEntry =
true;
786 if (characteristic != 0)
result =
result % characteristic;
787 s++;
m++; as++; am++;
819 const int columnIndex)
const
825 const int absoluteColumnIndex)
const
827 return getEntry(absoluteRowIndex, absoluteColumnIndex) ==
NULL;
834 string s =
"PolyMinorProcessor:";
839 int myIndexArray[500];
840 s +=
"\n considered submatrix has row indices: ";
844 if (
k != 0)
s +=
", ";
845 sprintf(
h,
"%d", myIndexArray[
k]);
s +=
h;
847 s +=
" (first row of matrix has index 0)";
848 s +=
"\n considered submatrix has column indices: ";
852 if (
k != 0)
s +=
", ";
853 sprintf(
h,
"%d", myIndexArray[
k]);
s +=
h;
855 s +=
" (first column of matrix has index 0)";
856 s +=
"\n size of considered minor(s): ";
867 for (
int i = 0;
i < n;
i++)
873 const int numberOfColumns,
874 const poly* polyMatrix)
878 for (
int i = 0;
i < n;
i++)
882 _rows = numberOfRows;
891 for (
int i = 0;
i < n;
i++)
896 const int* rowIndices,
897 const int* columnIndices,
909 const int* rowIndices,
910 const int* columnIndices,
911 const char* algorithm,
917 if (strcmp(algorithm,
"Laplace") == 0)
919 else if (strcmp(algorithm,
"Bareiss") == 0)
933 if (strcmp(algorithm,
"Laplace") == 0)
935 else if (strcmp(algorithm,
"Bareiss") == 0)
977 int s = 0;
int m = 0;
int as = 0;
int am = 0;
981 bool hadNonZeroEntry =
false;
989 for (
int c = 0; c <
k; c++)
995 hadNonZeroEntry =
true;
1009#ifdef COUNT_AND_PRINT_OPERATIONS
1014 s++;
m++; as++, am++;
1029 for (
int r = 0; r <
k; r++)
1035 hadNonZeroEntry =
true;
1049#ifdef COUNT_AND_PRINT_OPERATIONS
1054 s++;
m++; as++, am++;
1060 if (hadNonZeroEntry)
1063#ifdef COUNT_AND_PRINT_OPERATIONS
1089 const bool multipleMinors,
1101 0, 0, 0, 0, -1, -1);
1111 int s = 0;
int m = 0;
int as = 0;
int am = 0;
1115 bool hadNonZeroEntry =
false;
1124 for (
int c = 0; c <
k; c++)
1130 hadNonZeroEntry =
true;
1163#ifdef COUNT_AND_PRINT_OPERATIONS
1168 s++;
m++; as++; am++;
1183 for (
int r = 0; r <
k; r++)
1189 hadNonZeroEntry =
true;
1221#ifdef COUNT_AND_PRINT_OPERATIONS
1226 s++;
m++; as++; am++;
1238 if (hadNonZeroEntry)
1241#ifdef COUNT_AND_PRINT_OPERATIONS
1267 poly a = f1; poly
b = f2;
1271 b = f1; a = f2; bLen = aLen;
1291#ifdef COUNT_AND_PRINT_OPERATIONS
1325 number &c5,
int p5Len)
1327#ifdef COUNT_AND_PRINT_OPERATIONS
1354 while (bucketLm !=
NULL)
1365#ifdef COUNT_AND_PRINT_OPERATIONS
1367 multsMonForDiv += p5Len;
1379 helperPoly = bucketLm;
1380 helperPoly->next = p1;
1397 int *theRows=(
int*)
omAlloc(
k*
sizeof(
int));
1399 int *theColumns=(
int*)
omAlloc(
k*
sizeof(
int));
1404 0, 0, 0, 0, -1, -1);
1412 poly* tempMatrix = (poly*)
omAlloc(
k *
k *
sizeof(poly));
1415 for (
int r = 0; r <
k; r++)
1416 for (
int c = 0; c <
k; c++)
1422 int *rowPermutation=(
int*)
omAlloc(
k*
sizeof(
int));
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++)
1434 int minComplexity = -1;
int complexity = 0;
int bestRow = -1;
1436 for (
int i = r;
i <
k;
i++)
1438 pp = tempMatrix[rowPermutation[
i] *
k + r];
1441 if (minComplexity == -1)
1449 while ((
pp !=
NULL) && (complexity < minComplexity))
1453 if (complexity < minComplexity)
1455 minComplexity = complexity;
1459 if (minComplexity <= 1)
break;
1468 pNormalize(tempMatrix[rowPermutation[bestRow] *
k + r]);
1472 int j = rowPermutation[bestRow];
1473 rowPermutation[bestRow] = rowPermutation[r];
1474 rowPermutation[r] =
j;
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++)
1485 for (
int v = 0;
v <
k;
v++)
1487 if ((
v < r) && (u >
v))
1491 w = tempMatrix[rowPermutation[u] *
k +
v];
1502 divisor = tempMatrix[rowPermutation[r - 1] *
k + r - 1];
1504 divisorLength =
pLength(divisor);
1507 for (
int rr = r + 1; rr <
k; rr++)
1508 for (
int cc = r + 1; cc <
k; cc++)
1512 tempMatrix[rowPermutation[r] *
k + r],
1513 tempMatrix[rowPermutation[r] *
k + cc],
1514 tempMatrix[rowPermutation[rr] *
k + r]);
1517 tempMatrix[rowPermutation[r] *
k + r],
1518 tempMatrix[rowPermutation[r] *
k + cc],
1519 tempMatrix[rowPermutation[rr] *
k + r],
1520 divisor, divisorLC, divisorLength);
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++)
1528 for (
int v = 0;
v <
k;
v++)
1530 if ((
v <
k - 1) && (u >
v))
1534 w = tempMatrix[rowPermutation[u] *
k +
v];
1542 poly
result = tempMatrix[rowPermutation[
k - 1] *
k +
k - 1];
1543 tempMatrix[rowPermutation[
k - 1] *
k +
k - 1]=
NULL;
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)
Class Cache is a template-implementation of a cache with arbitrary classes for representing keys and ...
bool put(const KeyClass &key, const ValueClass &value)
Inserts the pair (key --> value) in the cache.
bool hasKey(const KeyClass &key) const
Checks whether the cache contains a pair (k --> v) such that k equals the given key.
ValueClass getValue(const KeyClass &key) const
Returns the value for a given key.
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...
int getResult() const
Accessor for the private field _result.
Class MinorKey can be used for representing keys in a cache for sub-determinantes; see class Cache.
void getAbsoluteColumnIndices(int *const target) const
A method for retrieving the 0-based indices of all columns encoded in this MinorKey.
void selectFirstRows(const int k, const MinorKey &mk)
This method redefines the set of rows represented by this MinorKey.
bool selectNextColumns(const int k, const MinorKey &mk)
This method redefines the set of columns represented by this MinorKey.
void reset()
A method for deleting all entries of _rowKey and _columnKey.
void selectFirstColumns(const int k, const MinorKey &mk)
This method redefines the set of columns represented by this MinorKey.
MinorKey getSubMinorKey(const int absoluteEraseRowIndex, const int absoluteEraseColumnIndex) const
A method for retrieving a sub-MinorKey resulting from omitting one row and one column of this MinorKe...
int getAbsoluteColumnIndex(const int i) const
A method for retrieving the (0-based) index of the i-th column in the set of columns encoded in this.
void set(const int lengthOfRowArray, const unsigned int *rowKey, const int lengthOfColumnArray, const unsigned int *columnKey)
A setter method for class MinorKey.
void getAbsoluteRowIndices(int *const target) const
A method for retrieving the 0-based indices of all rows encoded in this MinorKey.
int getRelativeRowIndex(const int i) const
A method for retrieving the (0-based) relative index of the i-th row in this MinorKey.
bool selectNextRows(const int k, const MinorKey &mk)
This method redefines the set of rows represented by this MinorKey.
int getRelativeColumnIndex(const int i) const
A method for retrieving the (0-based) relative index of the i-th column in this MinorKey.
int compare(const MinorKey &mk) const
A comparator for two instances of MinorKey.
int getAbsoluteRowIndex(const int i) const
A method for retrieving the (0-based) index of the i-th row in the set of rows encoded in this.
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...
int getAdditions() const
A method for accessing the additions performed while computing this minor.
int getAccumulatedAdditions() const
A method for accessing the additions performed while computing this minor, including all nested addit...
int getMultiplications() const
A method for accessing the multiplications performed while computing this minor.
void incrementRetrievals()
A method for incrementing the number of performed retrievals of this instance of MinorValue.
int getAccumulatedMultiplications() const
A method for accessing the multiplications performed while computing this minor, including all nested...
~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...
poly getResult() const
Accessor for the private field _result.
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 ....
const CanonicalForm int s
const Variable & v
< [in] a sqrfree bivariate poly
void kBucketClear(kBucket_pt bucket, poly *p, int *length)
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...
void kBucketDestroy(kBucket_pt *bucket_pt)
kBucket_pt kBucketCreate(const ring bucket_ring)
Creation/Destruction of buckets.
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)
const poly kBucketGetLm(kBucket_pt bucket)
poly kNF(ideal F, ideal Q, poly p, int syzComp, int lazyReduce)
static number & pGetCoeff(poly p)
return an alias to the leading coefficient of p assumes that p != NULL NOTE: not copy
#define omFreeSize(addr, size)
static int pLength(poly a)
static poly p_Add_q(poly p, poly q, const ring r)
static poly p_Mult_q(poly p, poly q, const ring r)
static void p_ExpVectorSub(poly p1, poly p2, const ring r)
static poly pReverse(poly p)
static void p_Delete(poly *p, const ring r)
static poly pp_Mult_qq(poly p, poly q, const ring r)
VAR ring currRing
Widely used global variable which specifies the current polynomial ring for Singular interpreter and ...
Compatibility layer for legacy polynomial operations (over currRing)
#define pSetCoeff(p, n)
deletes old coeff before setting the new one
#define pCopy(p)
return a copy of the poly
void PrintS(const char *s)