42 for (
i = 1;
i <= rows;
i++ )
43 for (
j = 1;
j <= rows;
j++ )
52 int i,
j, rows =
M.rows(), cols =
M.columns();
53 for (
i = 1;
i <= rows;
i++ )
54 for (
j = 1;
j <= cols;
j++ )
65 else if ( oldpivot.
isZero() )
67 else if (
level( oldpivot ) >
level( newpivot ) )
69 else if (
level( oldpivot ) <
level( newpivot ) )
72 return ( newpivot.
lc() < oldpivot.
lc() );
90 if (
M(
j,
i) != 0 )
break;
91 if (
j >
nrows )
return false;
94 pivotrecip = 1 /
M(
i,
i);
99 if ( rowpivot == 0 )
continue;
101 M(
j,
k) -=
M(
i,
k) * rowpivot;
115 int rows =
M.rows(), cols =
M.columns();
123 for (
i = 0;
i < rows;
i++ ) {
124 mm[
i] =
new int[cols];
135 DEBOUT( cerr,
"trying prime(" << pno <<
") = " );
141 for (
i = 1;
i <= rows;
i++ )
142 for (
j = 1;
j <= cols;
j++ )
145 ok =
solve( mm, rows, cols );
151 for (
i = 1;
i <= rows;
i++ )
152 for (
j = rows+1;
j <= cols;
j++ )
153 MM(
i,
j) = mm[
i-1][
j-1];
160 DEBOUT( cerr,
"trying prime(" << pno <<
") = " );
165 for (
i = 1;
i <= rows;
i++ )
166 for (
j = 1;
j <= cols;
j++ )
169 ok =
solve( mm, rows, cols );
175 for (
i = 1;
i <= rows;
i++ )
176 for (
j = rows+1;
j <= cols;
j++ )
189 for (
i = 1;
i <= rows;
i++ ) {
190 for (
j = rows+1;
j <= cols;
j++ )
191 if ( MM(
i,
j) > Qhalf )
207 for (
i = 0;
i < rows && ok;
i++ )
208 for (
j = 0;
j < rows && ok;
j++ )
226 ASSERT( rows <=
M.rows() && rows <=
M.columns() && rows > 0,
"undefined determinant" );
229 else if ( rows == 2 )
230 return M(1,1)*
M(2,2)-
M(2,1)*
M(1,2);
235 int n,
i, intdet,
p, pno;
236 for (
i = 0;
i < rows;
i++ )
238 mm[
i] =
new int[rows];
278 for (
i = 0;
i < rows;
i++ )
288 for (
i = 1;
i <= rows;
i++ ) {
290 for (
j =
i+1;
j <= rows;
j++ ) {
296 if (
pivot.isZero() )
303 for (
j =
i+1;
j <= rows;
j++ )
310 for (
k =
i+1;
k <= rows;
k++ )
316 for (
i = 1;
i <= rows;
i++ )
318 return pivot / divisor;
327 ASSERT( rows <=
M.rows() && rows <=
M.columns() && rows > 0,
"undefined determinant" );
330 else if ( rows == 2 )
331 return M(1,1)*
M(2,2)-
M(2,1)*
M(1,2);
336 int i,
p, pcount, pno, intdet;
340 for (
i = 0;
i < rows;
i++ ) {
341 mm[
i] =
new int[rows];
419 for (
i = 0;
i < rows;
i++ )
428 for (
i = 1;
i <= rows;
i++ ) {
430 for (
j =
i+1;
j <= rows;
j++ ) {
436 if (
pivot.isZero() )
442 for (
j =
i+1;
j <= rows;
j++ ) {
447 for (
k =
i+1;
k <= rows;
k++ )
453 for (
i = 1;
i <= rows;
i++ )
455 return pivot / divisor;
463 int rows =
M.rows(), cols =
M.columns();
466 for (
i = 1;
i <= rows;
i++ )
467 for (
j = 1;
j <= rows;
j++ )
469 DEBOUTLN( cerr,
"bound(matrix)^2 = " << sum );
471 for (
j = rows+1;
j <= cols;
j++ ) {
473 for (
i = 1;
i <= rows;
i++ )
475 if ( vsum > vmax ) vmax = vsum;
477 DEBOUTLN( cerr,
"bound(lhs)^2 = " << vmax );
479 DEBOUTLN( cerr,
"bound(overall)^2 = " << sum );
481 return sqrt( sum ) + 1;
490 for (
i = 1;
i <= rows;
i++ ) {
492 for (
j = 1;
j <= rows;
j++ )
508 int rowpivot, pivotrecip;
516 if ( extmat[
j][
i] != 0 )
break;
523 swap = extmat[
i]; extmat[
i] = extmat[
j]; extmat[
j] =
swap;
525 pivotrecip =
ff_inv( extmat[
i][
i] );
528 rowi[
j] =
ff_mul( pivotrecip, rowi[
j] );
532 if ( rowpivot == 0 )
continue;
541 for (
j = 0;
j <
i;
j++ ) {
544 if ( rowpivot == 0 )
continue;
550 DEBOUTLN( cerr,
"solve succeeded" );
559 int divisor, multiplier, rowii, rowji;
567 for (
i = 0;
i < n;
i++ ) {
569 for (
j =
i;
j < n;
j++ )
570 if ( extmat[
j][
i] != 0 )
break;
571 if (
j == n )
return 0;
573 multiplier =
ff_neg( multiplier );
574 swap = extmat[
i]; extmat[
i] = extmat[
j]; extmat[
j] =
swap;
578 for (
j =
i+1;
j < n;
j++ ) {
581 if ( rowji == 0 )
continue;
582 divisor =
ff_mul( divisor, rowii );
583 for (
k =
i;
k < n;
k++ )
588 for (
i = 0;
i < n;
i++ )
589 multiplier =
ff_mul( multiplier, extmat[
i][
i] );
600 for (
i = 1;
i <= n;
i++ )
602 for (
i = 1;
i <= n;
i++ ) {
603 q =
Q / ( z - a[
i] );
604 p = q / q( a[
i], z );
606 for (
j =
p;
j.hasTerms();
j++ )
607 x[
i] +=
w[
j.exp()+1] *
j.coeff();
declarations of higher level algorithms.
void FACTORY_PUBLIC chineseRemainder(const CanonicalForm &x1, const CanonicalForm &q1, const CanonicalForm &x2, const CanonicalForm &q2, CanonicalForm &xnew, CanonicalForm &qnew)
void chineseRemainder ( const CanonicalForm & x1, const CanonicalForm & q1, const CanonicalForm & x2,...
#define ASSERT(expression, message)
Iterators for CanonicalForm's.
bool matrix_in_Z(const CFMatrix &M, int rows)
bool betterpivot(const CanonicalForm &oldpivot, const CanonicalForm &newpivot)
CanonicalForm detbound(const CFMatrix &M, int rows)
void solveVandermondeT(const CFArray &a, const CFArray &w, CFArray &x, const Variable &z)
static CanonicalForm bound(const CFMatrix &M)
bool linearSystemSolve(CFMatrix &M)
int determinant(int **extmat, int n)
bool solve(int **extmat, int nrows, int ncols)
CanonicalForm determinant2(const CFMatrix &M, int rows)
static bool fill_int_mat(const CFMatrix &M, int **m, int rows)
int cf_getBigPrime(int i)
class to iterate through CanonicalForm's
factory's class for variables
functions to print debug output
#define DEBINCLEVEL(stream, msg)
#define DEBOUTENDL(stream)
#define DEBOUTLN(stream, objects)
#define DEBOUT(stream, objects)
#define DEBDECLEVEL(stream, msg)
bool isZero(const CFArray &A)
checks if entries of A are zero
operations in a finite prime field F_p.
int ff_mul(const int a, const int b)
int ff_sub(const int a, const int b)
bool pivot(const matrix aMat, const int r1, const int r2, const int c1, const int c2, int *bestR, int *bestC, const ring R)
This code computes a score for each non-zero matrix entry in aMat[r1..r2, c1..c2].
gmp_float sqrt(const gmp_float &a)
#define TIMING_DEFINE_PRINT(t)