Functions | Variables
cf_linsys.cc File Reference

solve linear systems and compute determinants of matrices More...

#include "config.h"
#include "cf_assert.h"
#include "debug.h"
#include "timing.h"
#include "cf_defs.h"
#include "cf_primes.h"
#include "canonicalform.h"
#include "cf_iter.h"
#include "cf_algorithm.h"
#include "ffops.h"

Go to the source code of this file.

Functions

 TIMING_DEFINE_PRINT (det_mapping) TIMING_DEFINE_PRINT(det_determinant) TIMING_DEFINE_PRINT(det_chinese) TIMING_DEFINE_PRINT(det_bound) TIMING_DEFINE_PRINT(det_numprimes) static bool solve(int **extmat
 
int determinant (int **extmat, int n)
 
static CanonicalForm bound (const CFMatrix &M)
 
CanonicalForm detbound (const CFMatrix &M, int rows)
 
bool matrix_in_Z (const CFMatrix &M, int rows)
 
bool matrix_in_Z (const CFMatrix &M)
 
bool betterpivot (const CanonicalForm &oldpivot, const CanonicalForm &newpivot)
 
bool linearSystemSolve (CFMatrix &M)
 
static bool fill_int_mat (const CFMatrix &M, int **m, int rows)
 
CanonicalForm determinant (const CFMatrix &M, int rows)
 
CanonicalForm determinant2 (const CFMatrix &M, int rows)
 
bool solve (int **extmat, int nrows, int ncols)
 
void solveVandermondeT (const CFArray &a, const CFArray &w, CFArray &x, const Variable &z)
 

Variables

int nrows
 
int int ncols
 
bool fuzzy_result
 

Detailed Description

solve linear systems and compute determinants of matrices

Definition in file cf_linsys.cc.

Function Documentation

§ betterpivot()

bool betterpivot ( const CanonicalForm oldpivot,
const CanonicalForm newpivot 
)

Definition at line 61 of file cf_linsys.cc.

62 {
63  if ( newpivot.isZero() )
64  return false;
65  else if ( oldpivot.isZero() )
66  return true;
67  else if ( level( oldpivot ) > level( newpivot ) )
68  return true;
69  else if ( level( oldpivot ) < level( newpivot ) )
70  return false;
71  else
72  return ( newpivot.lc() < oldpivot.lc() );
73 }
int level(const CanonicalForm &f)
CF_NO_INLINE bool isZero() const
Definition: cf_inline.cc:372
CanonicalForm lc() const
CanonicalForm CanonicalForm::lc (), Lc (), LC (), LC ( v ) const.

§ bound()

static CanonicalForm bound ( const CFMatrix M)
static

Definition at line 460 of file cf_linsys.cc.

461 {
462  DEBINCLEVEL( cerr, "bound" );
463  int rows = M.rows(), cols = M.columns();
464  CanonicalForm sum = 0;
465  int i, j;
466  for ( i = 1; i <= rows; i++ )
467  for ( j = 1; j <= rows; j++ )
468  sum += M(i,j) * M(i,j);
469  DEBOUTLN( cerr, "bound(matrix)^2 = " << sum );
470  CanonicalForm vmax = 0, vsum;
471  for ( j = rows+1; j <= cols; j++ ) {
472  vsum = 0;
473  for ( i = 1; i <= rows; i++ )
474  vsum += M(i,j) * M(i,j);
475  if ( vsum > vmax ) vmax = vsum;
476  }
477  DEBOUTLN( cerr, "bound(lhs)^2 = " << vmax );
478  sum += vmax;
479  DEBOUTLN( cerr, "bound(overall)^2 = " << sum );
480  DEBDECLEVEL( cerr, "bound" );
481  return sqrt( sum ) + 1;
482 }
#define DEBINCLEVEL(stream, msg)
Definition: debug.h:44
factory&#39;s main class
Definition: canonicalform.h:75
#define DEBDECLEVEL(stream, msg)
Definition: debug.h:45
#define M
Definition: sirandom.c:24
int j
Definition: myNF.cc:70
int rows() const
Definition: ftmpl_matrix.h:45
int columns() const
Definition: ftmpl_matrix.h:46
gmp_float sqrt(const gmp_float &a)
Definition: mpr_complex.cc:329
int i
Definition: cfEzgcd.cc:123
#define DEBOUTLN(stream, objects)
Definition: debug.h:49

§ detbound()

CanonicalForm detbound ( const CFMatrix M,
int  rows 
)

Definition at line 486 of file cf_linsys.cc.

487 {
488  CanonicalForm sum = 0, prod = 2;
489  int i, j;
490  for ( i = 1; i <= rows; i++ ) {
491  sum = 0;
492  for ( j = 1; j <= rows; j++ )
493  sum += M(i,j) * M(i,j);
494  prod *= 1 + sqrt(sum);
495  }
496  return prod;
497 }
factory&#39;s main class
Definition: canonicalform.h:75
#define M
Definition: sirandom.c:24
int j
Definition: myNF.cc:70
gmp_float sqrt(const gmp_float &a)
Definition: mpr_complex.cc:329
int i
Definition: cfEzgcd.cc:123
fq_nmod_poly_t prod
Definition: facHensel.cc:95

§ determinant() [1/2]

int determinant ( int **  extmat,
int  n 
)

Definition at line 556 of file cf_linsys.cc.

557 {
558  int i, j, k;
559  int divisor, multiplier, rowii, rowji; // all FF
560  int * rowi; // FF
561  int * rowj; // FF
562  int * swap; // FF
563  // triangularization
564  multiplier = 1;
565  divisor = 1;
566 
567  for ( i = 0; i < n; i++ ) {
568  //find "pivot"
569  for (j = i; j < n; j++ )
570  if ( extmat[j][i] != 0 ) break;
571  if ( j == n ) return 0;
572  if ( j != i ) {
573  multiplier = ff_neg( multiplier );
574  swap = extmat[i]; extmat[i] = extmat[j]; extmat[j] = swap;
575  }
576  rowi = extmat[i];
577  rowii = rowi[i];
578  for ( j = i+1; j < n; j++ ) {
579  rowj = extmat[j];
580  rowji = rowj[i];
581  if ( rowji == 0 ) continue;
582  divisor = ff_mul( divisor, rowii );
583  for ( k = i; k < n; k++ )
584  rowj[k] = ff_sub( ff_mul( rowj[k], rowii ), ff_mul( rowi[k], rowji ) );
585  }
586  }
587  multiplier = ff_mul( multiplier, ff_inv( divisor ) );
588  for ( i = 0; i < n; i++ )
589  multiplier = ff_mul( multiplier, extmat[i][i] );
590  return multiplier;
591 }
int ff_mul(const int a, const int b)
Definition: ffops.h:141
int k
Definition: cfEzgcd.cc:93
int j
Definition: myNF.cc:70
int ff_sub(const int a, const int b)
Definition: ffops.h:114
int ff_neg(const int a)
Definition: ffops.h:128
int i
Definition: cfEzgcd.cc:123
#define swap(_i, _j)
int ff_inv(const int a)
Definition: ffops.h:149

§ determinant() [2/2]

CanonicalForm determinant ( const CFMatrix M,
int  rows 
)

Definition at line 222 of file cf_linsys.cc.

223 {
224  typedef int* int_ptr;
225 
226  ASSERT( rows <= M.rows() && rows <= M.columns() && rows > 0, "undefined determinant" );
227  if ( rows == 1 )
228  return M(1,1);
229  else if ( rows == 2 )
230  return M(1,1)*M(2,2)-M(2,1)*M(1,2);
231  else if ( matrix_in_Z( M, rows ) )
232  {
233  int ** mm = new int_ptr[rows];
234  CanonicalForm x, q, Qhalf, B;
235  int n, i, intdet, p, pno;
236  for ( i = 0; i < rows; i++ )
237  {
238  mm[i] = new int[rows];
239  }
240  pno = 0; n = 0;
241  TIMING_START(det_bound);
242  B = detbound( M, rows );
243  TIMING_END(det_bound);
244  q = 1;
245  TIMING_START(det_numprimes);
246  while ( B > q && n < cf_getNumBigPrimes() )
247  {
248  q *= cf_getBigPrime( n );
249  n++;
250  }
251  TIMING_END(det_numprimes);
252 
253  CFArray X(1,n), Q(1,n);
254 
255  while ( pno < n )
256  {
257  p = cf_getBigPrime( pno );
258  setCharacteristic( p );
259  // map matrix into char p
260  TIMING_START(det_mapping);
261  fill_int_mat( M, mm, rows );
262  TIMING_END(det_mapping);
263  pno++;
264  DEBOUT( cerr, "." );
265  TIMING_START(det_determinant);
266  intdet = determinant( mm, rows );
267  TIMING_END(det_determinant);
268  setCharacteristic( 0 );
269  X[pno] = intdet;
270  Q[pno] = p;
271  }
272  TIMING_START(det_chinese);
273  chineseRemainder( X, Q, x, q );
274  TIMING_END(det_chinese);
275  Qhalf = q / 2;
276  if ( x > Qhalf )
277  x = x - q;
278  for ( i = 0; i < rows; i++ )
279  delete [] mm[i];
280  delete [] mm;
281  return x;
282  }
283  else
284  {
285  CFMatrix m( M );
286  CanonicalForm divisor = 1, pivot, mji;
287  int i, j, k, sign = 1;
288  for ( i = 1; i <= rows; i++ ) {
289  pivot = m(i,i); k = i;
290  for ( j = i+1; j <= rows; j++ ) {
291  if ( betterpivot( pivot, m(j,i) ) ) {
292  pivot = m(j,i);
293  k = j;
294  }
295  }
296  if ( pivot.isZero() )
297  return 0;
298  if ( i != k )
299  {
300  m.swapRow( i, k );
301  sign = -sign;
302  }
303  for ( j = i+1; j <= rows; j++ )
304  {
305  if ( ! m(j,i).isZero() )
306  {
307  divisor *= pivot;
308  mji = m(j,i);
309  m(j,i) = 0;
310  for ( k = i+1; k <= rows; k++ )
311  m(j,k) = m(j,k) * pivot - m(i,k)*mji;
312  }
313  }
314  }
315  pivot = sign;
316  for ( i = 1; i <= rows; i++ )
317  pivot *= m(i,i);
318  return pivot / divisor;
319  }
320 }
int cf_getNumBigPrimes()
Definition: cf_primes.cc:45
CanonicalForm detbound(const CFMatrix &M, int rows)
Definition: cf_linsys.cc:486
bool betterpivot(const CanonicalForm &oldpivot, const CanonicalForm &newpivot)
Definition: cf_linsys.cc:61
return P p
Definition: myNF.cc:203
#define DEBOUT(stream, objects)
Definition: debug.h:47
factory&#39;s main class
Definition: canonicalform.h:75
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].
int k
Definition: cfEzgcd.cc:93
#define Q
Definition: sirandom.c:25
int determinant(int **extmat, int n)
Definition: cf_linsys.cc:556
void setCharacteristic(int c)
Definition: cf_char.cc:23
#define M
Definition: sirandom.c:24
#define TIMING_END(t)
Definition: timing.h:88
int j
Definition: myNF.cc:70
static bool fill_int_mat(const CFMatrix &M, int **m, int rows)
Definition: cf_linsys.cc:203
int rows() const
Definition: ftmpl_matrix.h:45
int columns() const
Definition: ftmpl_matrix.h:46
int m
Definition: cfEzgcd.cc:119
int i
Definition: cfEzgcd.cc:123
void 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...
Definition: cf_chinese.cc:52
b *CanonicalForm B
Definition: facBivar.cc:51
#define TIMING_START(t)
Definition: timing.h:87
int cf_getBigPrime(int i)
Definition: cf_primes.cc:39
Variable x
Definition: cfModGcd.cc:4023
bool isZero(const CFArray &A)
checks if entries of A are zero
#define ASSERT(expression, message)
Definition: cf_assert.h:99
bool matrix_in_Z(const CFMatrix &M, int rows)
Definition: cf_linsys.cc:39
int * int_ptr
Definition: structs.h:57
static int sign(int x)
Definition: ring.cc:3328

§ determinant2()

CanonicalForm determinant2 ( const CFMatrix M,
int  rows 
)

Definition at line 323 of file cf_linsys.cc.

324 {
325  typedef int* int_ptr;
326 
327  ASSERT( rows <= M.rows() && rows <= M.columns() && rows > 0, "undefined determinant" );
328  if ( rows == 1 )
329  return M(1,1);
330  else if ( rows == 2 )
331  return M(1,1)*M(2,2)-M(2,1)*M(1,2);
332  else if ( matrix_in_Z( M, rows ) ) {
333  int ** mm = new int_ptr[rows];
334  CanonicalForm QQ, Q, Qhalf, mnew, q, qnew, B;
335  CanonicalForm det, detnew, qdet;
336  int i, p, pcount, pno, intdet;
337  bool ok;
338 
339  // initialize room to hold the result and the result mod p
340  for ( i = 0; i < rows; i++ ) {
341  mm[i] = new int[rows];
342  }
343 
344  // calculate the bound for the result
345  B = detbound( M, rows );
346 
347  // find a first solution mod p
348  pno = 0;
349  do {
350  p = cf_getBigPrime( pno );
351  setCharacteristic( p );
352  // map matrix into char p
353  ok = fill_int_mat( M, mm, rows );
354  pno++;
355  } while ( ! ok && pno < cf_getNumPrimes() );
356  // initialize the result matrix with first solution
357  // solve mod p
358  DEBOUT( cerr, "." );
359  intdet = determinant( mm, rows );
360  setCharacteristic( 0 );
361  det = intdet;
362  // Q so far
363  Q = p;
364  QQ = p;
365  while ( Q < B && cf_getNumPrimes() > pno ) {
366  // find a first solution mod p
367  do {
368  p = cf_getBigPrime( pno );
369  setCharacteristic( p );
370  // map matrix into char p
371  ok = fill_int_mat( M, mm, rows );
372  pno++;
373  } while ( ! ok && pno < cf_getNumPrimes() );
374  // initialize the result matrix with first solution
375  // solve mod p
376  DEBOUT( cerr, "." );
377  intdet = determinant( mm, rows );
378  setCharacteristic( 0 );
379  qdet = intdet;
380  // Q so far
381  q = p;
382  QQ *= p;
383  pcount = 0;
384  while ( QQ < B && cf_getNumPrimes() > pno && pcount < 500 ) {
385  do {
386  p = cf_getBigPrime( pno );
387  setCharacteristic( p );
388  ok = true;
389  // map matrix into char p
390  ok = fill_int_mat( M, mm, rows );
391  pno++;
392  } while ( ! ok && cf_getNumPrimes() > pno );
393  // solve mod p
394  DEBOUT( cerr, "." );
395  intdet = determinant( mm, rows );
396  // found a solution mod p
397  // now chinese remainder it to a solution mod Q*p
398  setCharacteristic( 0 );
399  chineseRemainder( qdet, q, intdet, p, detnew, qnew );
400  qdet = detnew;
401  q = qnew;
402  QQ *= p;
403  pcount++;
404  }
405  DEBOUT( cerr, "*" );
406  chineseRemainder( det, Q, qdet, q, detnew, qnew );
407  Q = qnew;
408  QQ = Q;
409  det = detnew;
410  }
411  if ( ! ok )
412  fuzzy_result = true;
413  else
414  fuzzy_result = false;
415  // store the result in M
416  Qhalf = Q / 2;
417  if ( det > Qhalf )
418  det = det - Q;
419  for ( i = 0; i < rows; i++ )
420  delete [] mm[i];
421  delete [] mm;
422  return det;
423  }
424  else {
425  CFMatrix m( M );
426  CanonicalForm divisor = 1, pivot, mji;
427  int i, j, k, sign = 1;
428  for ( i = 1; i <= rows; i++ ) {
429  pivot = m(i,i); k = i;
430  for ( j = i+1; j <= rows; j++ ) {
431  if ( betterpivot( pivot, m(j,i) ) ) {
432  pivot = m(j,i);
433  k = j;
434  }
435  }
436  if ( pivot.isZero() )
437  return 0;
438  if ( i != k ) {
439  m.swapRow( i, k );
440  sign = -sign;
441  }
442  for ( j = i+1; j <= rows; j++ ) {
443  if ( ! m(j,i).isZero() ) {
444  divisor *= pivot;
445  mji = m(j,i);
446  m(j,i) = 0;
447  for ( k = i+1; k <= rows; k++ )
448  m(j,k) = m(j,k) * pivot - m(i,k)*mji;
449  }
450  }
451  }
452  pivot = sign;
453  for ( i = 1; i <= rows; i++ )
454  pivot *= m(i,i);
455  return pivot / divisor;
456  }
457 }
CanonicalForm detbound(const CFMatrix &M, int rows)
Definition: cf_linsys.cc:486
bool betterpivot(const CanonicalForm &oldpivot, const CanonicalForm &newpivot)
Definition: cf_linsys.cc:61
return P p
Definition: myNF.cc:203
#define DEBOUT(stream, objects)
Definition: debug.h:47
bool fuzzy_result
Definition: cf_linsys.cc:75
factory&#39;s main class
Definition: canonicalform.h:75
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].
int k
Definition: cfEzgcd.cc:93
#define Q
Definition: sirandom.c:25
int determinant(int **extmat, int n)
Definition: cf_linsys.cc:556
void setCharacteristic(int c)
Definition: cf_char.cc:23
#define M
Definition: sirandom.c:24
int j
Definition: myNF.cc:70
static bool fill_int_mat(const CFMatrix &M, int **m, int rows)
Definition: cf_linsys.cc:203
int rows() const
Definition: ftmpl_matrix.h:45
int columns() const
Definition: ftmpl_matrix.h:46
int m
Definition: cfEzgcd.cc:119
int i
Definition: cfEzgcd.cc:123
void 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...
Definition: cf_chinese.cc:52
b *CanonicalForm B
Definition: facBivar.cc:51
int cf_getBigPrime(int i)
Definition: cf_primes.cc:39
bool isZero(const CFArray &A)
checks if entries of A are zero
#define ASSERT(expression, message)
Definition: cf_assert.h:99
bool matrix_in_Z(const CFMatrix &M, int rows)
Definition: cf_linsys.cc:39
int * int_ptr
Definition: structs.h:57
static int sign(int x)
Definition: ring.cc:3328
int cf_getNumPrimes()
Definition: cf_primes.cc:23

§ fill_int_mat()

static bool fill_int_mat ( const CFMatrix M,
int **  m,
int  rows 
)
static

Definition at line 203 of file cf_linsys.cc.

204 {
205  int i, j;
206  bool ok = true;
207  for ( i = 0; i < rows && ok; i++ )
208  for ( j = 0; j < rows && ok; j++ )
209  {
210  if ( M(i+1,j+1).isZero() )
211  m[i][j] = 0;
212  else
213  {
214  m[i][j] = mapinto( M(i+1,j+1) ).intval();
215 // ok = m[i][j] != 0;
216  }
217  }
218  return ok;
219 }
long intval() const
conversion functions
#define M
Definition: sirandom.c:24
int j
Definition: myNF.cc:70
int m
Definition: cfEzgcd.cc:119
int i
Definition: cfEzgcd.cc:123
CanonicalForm mapinto(const CanonicalForm &f)
bool isZero(const CFArray &A)
checks if entries of A are zero

§ linearSystemSolve()

bool linearSystemSolve ( CFMatrix M)

Definition at line 78 of file cf_linsys.cc.

79 {
80  typedef int* int_ptr;
81 
82  if ( ! matrix_in_Z( M ) ) {
83  int nrows = M.rows(), ncols = M.columns();
84  int i, j, k;
85  CanonicalForm rowpivot, pivotrecip;
86  // triangularization
87  for ( i = 1; i <= nrows; i++ ) {
88  //find "pivot"
89  for (j = i; j <= nrows; j++ )
90  if ( M(j,i) != 0 ) break;
91  if ( j > nrows ) return false;
92  if ( j != i )
93  M.swapRow( i, j );
94  pivotrecip = 1 / M(i,i);
95  for ( j = 1; j <= ncols; j++ )
96  M(i,j) *= pivotrecip;
97  for ( j = i+1; j <= nrows; j++ ) {
98  rowpivot = M(j,i);
99  if ( rowpivot == 0 ) continue;
100  for ( k = i; k <= ncols; k++ )
101  M(j,k) -= M(i,k) * rowpivot;
102  }
103  }
104  // matrix is now upper triangular with 1s down the diagonal
105  // back-substitute
106  for ( i = nrows-1; i > 0; i-- ) {
107  for ( j = nrows+1; j <= ncols; j++ ) {
108  for ( k = i+1; k <= nrows; k++ )
109  M(i,j) -= M(k,j) * M(i,k);
110  }
111  }
112  return true;
113  }
114  else {
115  int rows = M.rows(), cols = M.columns();
116  CFMatrix MM( rows, cols );
117  int ** mm = new int_ptr[rows];
118  CanonicalForm Q, Qhalf, mnew, qnew, B;
119  int i, j, p, pno;
120  bool ok;
121 
122  // initialize room to hold the result and the result mod p
123  for ( i = 0; i < rows; i++ ) {
124  mm[i] = new int[cols];
125  }
126 
127  // calculate the bound for the result
128  B = bound( M );
129  DEBOUTLN( cerr, "bound = " << B );
130 
131  // find a first solution mod p
132  pno = 0;
133  do {
134  DEBOUTSL( cerr );
135  DEBOUT( cerr, "trying prime(" << pno << ") = " );
136  p = cf_getBigPrime( pno );
137  DEBOUT( cerr, p );
138  DEBOUTENDL( cerr );
139  setCharacteristic( p );
140  // map matrix into char p
141  for ( i = 1; i <= rows; i++ )
142  for ( j = 1; j <= cols; j++ )
143  mm[i-1][j-1] = mapinto( M(i,j) ).intval();
144  // solve mod p
145  ok = solve( mm, rows, cols );
146  pno++;
147  } while ( ! ok );
148 
149  // initialize the result matrix with first solution
150  setCharacteristic( 0 );
151  for ( i = 1; i <= rows; i++ )
152  for ( j = rows+1; j <= cols; j++ )
153  MM(i,j) = mm[i-1][j-1];
154 
155  // Q so far
156  Q = p;
157  while ( Q < B && pno < cf_getNumBigPrimes() ) {
158  do {
159  DEBOUTSL( cerr );
160  DEBOUT( cerr, "trying prime(" << pno << ") = " );
161  p = cf_getBigPrime( pno );
162  DEBOUT( cerr, p );
163  DEBOUTENDL( cerr );
164  setCharacteristic( p );
165  for ( i = 1; i <= rows; i++ )
166  for ( j = 1; j <= cols; j++ )
167  mm[i-1][j-1] = mapinto( M(i,j) ).intval();
168  // solve mod p
169  ok = solve( mm, rows, cols );
170  pno++;
171  } while ( ! ok );
172  // found a solution mod p
173  // now chinese remainder it to a solution mod Q*p
174  setCharacteristic( 0 );
175  for ( i = 1; i <= rows; i++ )
176  for ( j = rows+1; j <= cols; j++ )
177  {
178  chineseRemainder( MM[i][j], Q, CanonicalForm(mm[i-1][j-1]), CanonicalForm(p), mnew, qnew );
179  MM(i, j) = mnew;
180  }
181  Q = qnew;
182  }
183  if ( pno == cf_getNumBigPrimes() )
184  fuzzy_result = true;
185  else
186  fuzzy_result = false;
187  // store the result in M
188  Qhalf = Q / 2;
189  for ( i = 1; i <= rows; i++ ) {
190  for ( j = rows+1; j <= cols; j++ )
191  if ( MM(i,j) > Qhalf )
192  M(i,j) = MM(i,j) - Q;
193  else
194  M(i,j) = MM(i,j);
195  delete [] mm[i-1];
196  }
197  delete [] mm;
198  return ! fuzzy_result;
199  }
200 }
long intval() const
conversion functions
int cf_getNumBigPrimes()
Definition: cf_primes.cc:45
static CanonicalForm bound(const CFMatrix &M)
Definition: cf_linsys.cc:460
return P p
Definition: myNF.cc:203
#define DEBOUT(stream, objects)
Definition: debug.h:47
bool fuzzy_result
Definition: cf_linsys.cc:75
#define DEBOUTSL(stream)
Definition: debug.h:46
factory&#39;s main class
Definition: canonicalform.h:75
bool solve(int **extmat, int nrows, int ncols)
Definition: cf_linsys.cc:504
int k
Definition: cfEzgcd.cc:93
#define Q
Definition: sirandom.c:25
void setCharacteristic(int c)
Definition: cf_char.cc:23
#define M
Definition: sirandom.c:24
#define DEBOUTENDL(stream)
Definition: debug.h:48
int j
Definition: myNF.cc:70
int nrows
Definition: cf_linsys.cc:32
int rows() const
Definition: ftmpl_matrix.h:45
int columns() const
Definition: ftmpl_matrix.h:46
int i
Definition: cfEzgcd.cc:123
CanonicalForm mapinto(const CanonicalForm &f)
void swapRow(int i, int j)
void 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...
Definition: cf_chinese.cc:52
b *CanonicalForm B
Definition: facBivar.cc:51
int cf_getBigPrime(int i)
Definition: cf_primes.cc:39
int int ncols
Definition: cf_linsys.cc:32
#define DEBOUTLN(stream, objects)
Definition: debug.h:49
bool matrix_in_Z(const CFMatrix &M, int rows)
Definition: cf_linsys.cc:39
int * int_ptr
Definition: structs.h:57

§ matrix_in_Z() [1/2]

bool matrix_in_Z ( const CFMatrix M,
int  rows 
)

Definition at line 39 of file cf_linsys.cc.

40 {
41  int i, j;
42  for ( i = 1; i <= rows; i++ )
43  for ( j = 1; j <= rows; j++ )
44  if ( ! M(i,j).inZ() )
45  return false;
46  return true;
47 }
#define M
Definition: sirandom.c:24
int j
Definition: myNF.cc:70
int i
Definition: cfEzgcd.cc:123

§ matrix_in_Z() [2/2]

bool matrix_in_Z ( const CFMatrix M)

Definition at line 50 of file cf_linsys.cc.

51 {
52  int i, j, rows = M.rows(), cols = M.columns();
53  for ( i = 1; i <= rows; i++ )
54  for ( j = 1; j <= cols; j++ )
55  if ( ! M(i,j).inZ() )
56  return false;
57  return true;
58 }
#define M
Definition: sirandom.c:24
int j
Definition: myNF.cc:70
int rows() const
Definition: ftmpl_matrix.h:45
int columns() const
Definition: ftmpl_matrix.h:46
int i
Definition: cfEzgcd.cc:123

§ solve()

bool solve ( int **  extmat,
int  nrows,
int  ncols 
)

Definition at line 504 of file cf_linsys.cc.

505 {
506  DEBINCLEVEL( cerr, "solve" );
507  int i, j, k;
508  int rowpivot, pivotrecip; // all FF
509  int * rowi; // FF
510  int * rowj; // FF
511  int * swap; // FF
512  // triangularization
513  for ( i = 0; i < nrows; i++ ) {
514  //find "pivot"
515  for (j = i; j < nrows; j++ )
516  if ( extmat[j][i] != 0 ) break;
517  if ( j == nrows ) {
518  DEBOUTLN( cerr, "solve failed" );
519  DEBDECLEVEL( cerr, "solve" );
520  return false;
521  }
522  if ( j != i ) {
523  swap = extmat[i]; extmat[i] = extmat[j]; extmat[j] = swap;
524  }
525  pivotrecip = ff_inv( extmat[i][i] );
526  rowi = extmat[i];
527  for ( j = 0; j < ncols; j++ )
528  rowi[j] = ff_mul( pivotrecip, rowi[j] );
529  for ( j = i+1; j < nrows; j++ ) {
530  rowj = extmat[j];
531  rowpivot = rowj[i];
532  if ( rowpivot == 0 ) continue;
533  for ( k = i; k < ncols; k++ )
534  rowj[k] = ff_sub( rowj[k], ff_mul( rowpivot, rowi[k] ) );
535  }
536  }
537  // matrix is now upper triangular with 1s down the diagonal
538  // back-substitute
539  for ( i = nrows-1; i >= 0; i-- ) {
540  rowi = extmat[i];
541  for ( j = 0; j < i; j++ ) {
542  rowj = extmat[j];
543  rowpivot = rowj[i];
544  if ( rowpivot == 0 ) continue;
545  for ( k = i; k < ncols; k++ )
546  rowj[k] = ff_sub( rowj[k], ff_mul( rowpivot, rowi[k] ) );
547  // for (k=nrows; k<ncols; k++) rowj[k] = ff_sub(rowj[k], ff_mul(rowpivot, rowi[k]));
548  }
549  }
550  DEBOUTLN( cerr, "solve succeeded" );
551  DEBDECLEVEL( cerr, "solve" );
552  return true;
553 }
#define DEBINCLEVEL(stream, msg)
Definition: debug.h:44
int ff_mul(const int a, const int b)
Definition: ffops.h:141
int k
Definition: cfEzgcd.cc:93
#define DEBDECLEVEL(stream, msg)
Definition: debug.h:45
int j
Definition: myNF.cc:70
int nrows
Definition: cf_linsys.cc:32
int ff_sub(const int a, const int b)
Definition: ffops.h:114
int i
Definition: cfEzgcd.cc:123
#define swap(_i, _j)
int int ncols
Definition: cf_linsys.cc:32
#define DEBOUTLN(stream, objects)
Definition: debug.h:49
int ff_inv(const int a)
Definition: ffops.h:149

§ solveVandermondeT()

void solveVandermondeT ( const CFArray a,
const CFArray w,
CFArray x,
const Variable z 
)

Definition at line 594 of file cf_linsys.cc.

595 {
596  CanonicalForm Q = 1, q, p;
597  CFIterator j;
598  int i, n = a.size();
599 
600  for ( i = 1; i <= n; i++ )
601  Q *= ( z - a[i] );
602  for ( i = 1; i <= n; i++ ) {
603  q = Q / ( z - a[i] );
604  p = q / q( a[i], z );
605  x[i] = 0;
606  for ( j = p; j.hasTerms(); j++ )
607  x[i] += w[j.exp()+1] * j.coeff();
608  }
609 }
return P p
Definition: myNF.cc:203
CF_NO_INLINE CanonicalForm coeff() const
get the current coefficient
factory&#39;s main class
Definition: canonicalform.h:75
#define Q
Definition: sirandom.c:25
int size() const
Definition: ftmpl_array.cc:92
int j
Definition: myNF.cc:70
int i
Definition: cfEzgcd.cc:123
class to iterate through CanonicalForm&#39;s
Definition: cf_iter.h:44
CF_NO_INLINE int hasTerms() const
check if iterator has reached < the end of CanonicalForm
CF_NO_INLINE int exp() const
get the current exponent

§ TIMING_DEFINE_PRINT()

TIMING_DEFINE_PRINT ( det_mapping  )

Variable Documentation

§ fuzzy_result

bool fuzzy_result

Definition at line 75 of file cf_linsys.cc.

§ ncols

int int ncols

Definition at line 32 of file cf_linsys.cc.

§ nrows

int nrows

Definition at line 32 of file cf_linsys.cc.