My Project
Macros | Functions
intvec.cc File Reference
#include "misc/auxiliary.h"
#include "misc/intvec.h"
#include "misc/options.h"
#include "omalloc/omalloc.h"

Go to the source code of this file.

Macros

#define INTVEC_CC
 

Functions

intvecivAdd (intvec *a, intvec *b)
 
intvecivAddShift (intvec *a, intvec *b, int s)
 
intvecivSub (intvec *a, intvec *b)
 
intvecivTranp (intvec *o)
 
int ivTrace (intvec *o)
 
intvecivMult (intvec *a, intvec *b)
 
static int ivColPivot (intvec *, int, int, int, int)
 
static void ivNegRow (intvec *, int)
 
static void ivSaveRow (intvec *, int)
 
static void ivSetRow (intvec *, int, int)
 
static void ivFreeRow (intvec *, int, int)
 
static void ivReduce (intvec *, int, int, int, int)
 
static void ivZeroElim (intvec *, int, int, int &)
 
static void ivRowContent (intvec *, int, int)
 
static void ivKernFromRow (intvec *, intvec *, intvec *, int, int, int)
 
static intvecivOptimizeKern (intvec *)
 
static int ivGcd (int, int)
 
static void ivOptRecursive (intvec *, intvec *, intvec *, int &, int &, int)
 
static void ivOptSolve (intvec *, intvec *, int &, int &)
 
static void ivContent (intvec *)
 
static int ivL1Norm (intvec *)
 
static int ivCondNumber (intvec *, int)
 
void ivTriangIntern (intvec *imat, int &ready, int &all)
 
intvecivSolveKern (intvec *imat, int dimtr)
 
intvecivConcat (intvec *a, intvec *b)
 

Macro Definition Documentation

◆ INTVEC_CC

#define INTVEC_CC

Definition at line 8 of file intvec.cc.

Function Documentation

◆ ivAdd()

intvec* ivAdd ( intvec a,
intvec b 
)

Definition at line 249 of file intvec.cc.

250 {
251  intvec * iv;
252  int mn, ma, i;
253  if (a->cols() != b->cols()) return NULL;
254  mn = si_min(a->rows(),b->rows());
255  ma = si_max(a->rows(),b->rows());
256  if (a->cols() == 1)
257  {
258  iv = new intvec(ma);
259  for (i=0; i<mn; i++) (*iv)[i] = (*a)[i] + (*b)[i];
260  if (ma > mn)
261  {
262  if (ma == a->rows())
263  {
264  for(i=mn; i<ma; i++) (*iv)[i] = (*a)[i];
265  }
266  else
267  {
268  for(i=mn; i<ma; i++) (*iv)[i] = (*b)[i];
269  }
270  }
271  return iv;
272  }
273  if (mn != ma) return NULL;
274  iv = new intvec(a);
275  for (i=0; i<mn*a->cols(); i++) { (*iv)[i] += (*b)[i]; }
276  return iv;
277 }
static int si_max(const int a, const int b)
Definition: auxiliary.h:124
static int si_min(const int a, const int b)
Definition: auxiliary.h:125
int i
Definition: cfEzgcd.cc:132
CanonicalForm b
Definition: cfModGcd.cc:4103
Definition: intvec.h:23
int cols() const
Definition: intvec.h:95
int rows() const
Definition: intvec.h:96
#define NULL
Definition: omList.c:12

◆ ivAddShift()

intvec* ivAddShift ( intvec a,
intvec b,
int  s 
)

Definition at line 279 of file intvec.cc.

280 {
281  if (a->cols() != b->cols()) return NULL;
282  if (a->cols() == 1)
283  {
284  int i;
285  intvec * iv;
286  int ma = si_max(a->rows(),b->rows()+s);
287  iv=new intvec(ma);
288  for (i=a->rows()-1;i>=0;i--)
289  (*iv)[i] = (*a)[i];
290  for(i=b->rows()+s-1;i>=s;i--)
291  (*iv)[i] += (*b)[i-s];
292  return iv;
293  }
294  return NULL;
295 }
const CanonicalForm int s
Definition: facAbsFact.cc:51

◆ ivColPivot()

static int ivColPivot ( intvec imat,
int  colpos,
int  rowpos,
int  ready,
int  all 
)
static

Definition at line 484 of file intvec.cc.

485 {
486  int rpiv;
487 
488  if (IMATELEM(*imat,rowpos,colpos)!=0)
489  return rowpos;
490  for (rpiv=ready+1;rpiv<=all;rpiv++)
491  {
492  if (IMATELEM(*imat,rpiv,colpos)!=0)
493  return rpiv;
494  }
495  return 0;
496 }
#define IMATELEM(M, I, J)
Definition: intvec.h:85

◆ ivConcat()

intvec* ivConcat ( intvec a,
intvec b 
)

Definition at line 822 of file intvec.cc.

823 {
824  int ac=a->cols();
825  int c = ac + b->cols(); int r = si_max(a->rows(),b->rows());
826  intvec * ab = new intvec(r,c,0);
827 
828  int i,j;
829  for (i=1; i<=a->rows(); i++)
830  {
831  for(j=1; j<=ac; j++)
832  IMATELEM(*ab,i,j) = IMATELEM(*a,i,j);
833  }
834  for (i=1; i<=b->rows(); i++)
835  {
836  for(j=1; j<=b->cols(); j++)
837  IMATELEM(*ab,i,j+ac) = IMATELEM(*b,i,j);
838  }
839  return ab;
840 }
int j
Definition: facHensel.cc:110

◆ ivCondNumber()

static int ivCondNumber ( intvec w,
int  l 
)
static

Definition at line 764 of file intvec.cc.

765 {
766  int l0=0, i;
767 
768  if (l<0)
769  {
770  for (i=w->rows()-1;i>=0;i--)
771  {
772  if ((*w)[i]<0) l0--;
773  }
774  if (l0==0)
775  {
776  for (i=w->rows()-1;i>=0;i--)
777  {
778  if ((*w)[i]>0) l0++;
779  }
780  }
781  return l0;
782  }
783  else
784  {
785  for (i=w->rows()-1;i>=0;i--)
786  {
787  if ((*w)[i]<0) return -1;
788  }
789  for (i=w->rows()-1;i>=0;i--)
790  {
791  if ((*w)[i]>0) l0++;
792  }
793  return l0;
794  }
795 }
int l
Definition: cfEzgcd.cc:100
const CanonicalForm & w
Definition: facAbsFact.cc:51

◆ ivContent()

static void ivContent ( intvec w)
static

Definition at line 797 of file intvec.cc.

798 {
799  int tgcd, m;
800  int i=w->rows()-1;
801 
802  loop
803  {
804  tgcd = (*w)[i--];
805  if (tgcd!=0) break;
806  if (i<0) return;
807  }
808  if (tgcd<0) tgcd = -tgcd;
809  if (tgcd==1) return;
810  loop
811  {
812  m = (*w)[i--];
813  if (m!=0) tgcd= ivGcd(tgcd, m);
814  if (tgcd==1) return;
815  if (i<0) break;
816  }
817  for (i=w->rows()-1;i>=0;i--)
818  (*w)[i] /= tgcd;
819 }
int m
Definition: cfEzgcd.cc:128
static int ivGcd(int, int)
Definition: intvec.cc:648
#define loop
Definition: structs.h:75

◆ ivFreeRow()

static void ivFreeRow ( intvec imat,
int  rowpos,
int  rpiv 
)
static

Definition at line 522 of file intvec.cc.

523 {
524  int i, j;
525 
526  for (j=rpiv-1;j>=rowpos;j--)
527  {
528  for (i=imat->cols();i!=0;i--)
529  IMATELEM(*imat,j+1,i) = IMATELEM(*imat,j,i);
530  }
531 }

◆ ivGcd()

static int ivGcd ( int  a,
int  b 
)
static

Definition at line 648 of file intvec.cc.

649 {
650  int x;
651 
652  if (a<0) a=-a;
653  if (b<0) b=-b;
654  if (b>a)
655  {
656  x=b;
657  b=a;
658  a=x;
659  }
660  while (b!=0)
661  {
662  x = a % b;
663  a = b;
664  b = x;
665  }
666  return a;
667 }
Variable x
Definition: cfModGcd.cc:4082

◆ ivKernFromRow()

static void ivKernFromRow ( intvec kern,
intvec imat,
intvec perm,
int  pos,
int  r,
int  c 
)
static

Definition at line 613 of file intvec.cc.

615 {
616  int piv, cp, g, i, j, k, s;
617 
618  for (i=c;i>(*perm)[r];i--)
619  {
620  IMATELEM(*kern,pos,i) = 1;
621  for (j=r;j!=0;j--)
622  {
623  cp = (*perm)[j];
624  s=0;
625  for(k=c;k>cp;k--)
626  s += IMATELEM(*imat,j,k)*IMATELEM(*kern,pos,k);
627  if (s!=0)
628  {
629  piv = IMATELEM(*imat,j,cp);
630  g = ivGcd(piv,s);
631  if (g!=1)
632  {
633  s /= g;
634  piv /= g;
635  }
636  for(k=c;k>cp;k--)
637  IMATELEM(*kern,pos,k) *= piv;
638  IMATELEM(*kern,pos,cp) = -s;
639  ivRowContent(kern,pos,cp);
640  }
641  }
642  if (IMATELEM(*kern,pos,i)<0)
643  ivNegRow(kern,pos);
644  pos--;
645  }
646 }
int k
Definition: cfEzgcd.cc:99
g
Definition: cfModGcd.cc:4090
static void ivRowContent(intvec *, int, int)
Definition: intvec.cc:589
static void ivNegRow(intvec *, int)
Definition: intvec.cc:498

◆ ivL1Norm()

static int ivL1Norm ( intvec w)
static

Definition at line 749 of file intvec.cc.

750 {
751  int i, j, s = 0;
752 
753  for (i=w->rows()-1;i>=0;i--)
754  {
755  j = (*w)[i];
756  if (j>0)
757  s += j;
758  else
759  s -= j;
760  }
761  return s;
762 }

◆ ivMult()

intvec* ivMult ( intvec a,
intvec b 
)

Definition at line 349 of file intvec.cc.

350 {
351  int i, j, k, sum,
352  ra = a->rows(), ca = a->cols(),
353  rb = b->rows(), cb = b->cols();
354  intvec * iv;
355  if (ca != rb) return NULL;
356  iv = new intvec(ra, cb, 0);
357  for (i=0; i<ra; i++)
358  {
359  for (j=0; j<cb; j++)
360  {
361  sum = 0;
362  for (k=0; k<ca; k++)
363  sum += (*a)[i*ca+k]*(*b)[k*cb+j];
364  (*iv)[i*cb+j] = sum;
365  }
366  }
367  return iv;
368 }

◆ ivNegRow()

static void ivNegRow ( intvec imat,
int  rpiv 
)
static

Definition at line 498 of file intvec.cc.

499 {
500  int i;
501  for (i=imat->cols();i!=0;i--)
502  IMATELEM(*imat,rpiv,i) = -(IMATELEM(*imat,rpiv,i));
503 }

◆ ivOptimizeKern()

static intvec * ivOptimizeKern ( intvec kern)
static

Definition at line 669 of file intvec.cc.

670 {
671  int i,l,j,c=kern->cols(),r=kern->rows();
672  intvec *res=new intvec(c);
673 
674  if (TEST_OPT_PROT)
675  Warn(" %d linear independent solutions\n",r);
676  for (i=r;i>1;i--)
677  {
678  for (j=c;j>0;j--)
679  {
680  (*res)[j-1] += IMATELEM(*kern,i,j);
681  }
682  }
683  ivContent(res);
684  if (r<11)
685  {
686  l = ivCondNumber(res,-c);
687  j = ivL1Norm(res);
688  ivOptRecursive(res, NULL, kern, l, j, r);
689  }
690  return res;
691 }
#define Warn
Definition: emacs.cc:77
CanonicalForm res
Definition: facAbsFact.cc:60
static void ivOptRecursive(intvec *, intvec *, intvec *, int &, int &, int)
Definition: intvec.cc:693
static int ivL1Norm(intvec *)
Definition: intvec.cc:749
static void ivContent(intvec *)
Definition: intvec.cc:797
static int ivCondNumber(intvec *, int)
Definition: intvec.cc:764
#define TEST_OPT_PROT
Definition: options.h:103

◆ ivOptRecursive()

static void ivOptRecursive ( intvec res,
intvec w,
intvec kern,
int &  l,
int &  j,
int  pos 
)
static

Definition at line 693 of file intvec.cc.

695 {
696  int m, k, d;
697  intvec *h;
698 
699  d=kern->rows();
700  d=96/(d*d);
701  if (d<3) d=3;
702  if (w!=0)
703  h = new intvec(w);
704  else
705  h = new intvec(res->rows());
706  for (m=d;m>0;m--)
707  {
708  for(k=h->rows()-1;k>=0;k--)
709  (*h)[k] += IMATELEM(*kern,pos,k+1);
710  if(pos>1)
711  ivOptRecursive(res, h, kern, l, j, pos-1);
712  else
713  ivOptSolve(res, h, l, j);
714  }
715  delete h;
716  if (pos>1)
717  ivOptRecursive(res, w, kern, l, j, pos-1);
718  else if (w!=NULL)
719  ivOptSolve(res, w, l, j);
720 }
static void ivOptSolve(intvec *, intvec *, int &, int &)
Definition: intvec.cc:722
STATIC_VAR Poly * h
Definition: janet.cc:971

◆ ivOptSolve()

static void ivOptSolve ( intvec res,
intvec w,
int &  l,
int &  j 
)
static

Definition at line 722 of file intvec.cc.

723 {
724  int l0, j0, k;
725 
726  l0 = ivCondNumber(w, l);
727  if (l0==l)
728  {
729  ivContent(w);
730  j0 = ivL1Norm(w);
731  if(j0<j)
732  {
733  j = j0;
734  for(k=w->rows()-1;k>=0;k--)
735  (*res)[k] = (*w)[k];
736  }
737  return;
738  }
739  if(l0>l)
740  {
741  l = l0;
742  ivContent(w);
743  j = ivL1Norm(w);
744  for(k=w->rows()-1;k>=0;k--)
745  (*res)[k] = (*w)[k];
746  }
747 }

◆ ivReduce()

static void ivReduce ( intvec imat,
int  rpiv,
int  colpos,
int  ready,
int  all 
)
static

Definition at line 533 of file intvec.cc.

535 {
536  int tgcd, ce, m1, m2, j, i;
537  int piv = IMATELEM(*imat,rpiv,colpos);
538 
539  for (j=all;j>ready;j--)
540  {
541  ivRowContent(imat, j, 1);
542  ce = IMATELEM(*imat,j,colpos);
543  if (ce!=0)
544  {
545  IMATELEM(*imat,j,colpos) = 0;
546  m1 = piv;
547  m2 = ce;
548  tgcd = ivGcd(m1, m2);
549  if (tgcd != 1)
550  {
551  m1 /= tgcd;
552  m2 /= tgcd;
553  }
554  for (i=imat->cols();i>colpos;i--)
555  {
556  IMATELEM(*imat,j,i) = IMATELEM(*imat,j,i)*m1-
557  IMATELEM(*imat,rpiv,i)*m2;
558  }
559  ivRowContent(imat, j, colpos+1);
560  }
561  }
562 }

◆ ivRowContent()

static void ivRowContent ( intvec imat,
int  rowpos,
int  colpos 
)
static

Definition at line 589 of file intvec.cc.

590 {
591  int tgcd, m;
592  int i=imat->cols();
593 
594  loop
595  {
596  tgcd = IMATELEM(*imat,rowpos,i--);
597  if (tgcd!=0) break;
598  if (i<colpos) return;
599  }
600  if (tgcd<0) tgcd = -tgcd;
601  if (tgcd==1) return;
602  loop
603  {
604  m = IMATELEM(*imat,rowpos,i--);
605  if (m!=0) tgcd= ivGcd(tgcd, m);
606  if (tgcd==1) return;
607  if (i<colpos) break;
608  }
609  for (i=imat->cols();i>=colpos;i--)
610  IMATELEM(*imat,rowpos,i) /= tgcd;
611 }

◆ ivSaveRow()

static void ivSaveRow ( intvec imat,
int  rpiv 
)
static

Definition at line 505 of file intvec.cc.

506 {
507  int i, j=imat->rows();
508 
509  for (i=imat->cols();i!=0;i--)
510  IMATELEM(*imat,j,i) = IMATELEM(*imat,rpiv,i);
511 }

◆ ivSetRow()

static void ivSetRow ( intvec imat,
int  rowpos,
int  colpos 
)
static

Definition at line 513 of file intvec.cc.

514 {
515  int i, j=imat->rows();
516 
517  for (i=imat->cols();i!=0;i--)
518  IMATELEM(*imat,rowpos,i) = IMATELEM(*imat,j,i);
519  ivRowContent(imat, rowpos, colpos);
520 }

◆ ivSolveKern()

intvec* ivSolveKern ( intvec imat,
int  dimtr 
)

Definition at line 442 of file intvec.cc.

443 {
444  int d=imat->cols();
445  int kdim=d-dimtr;
446  intvec *perm = new intvec(dimtr+1);
447  intvec *kern = new intvec(kdim,d,0);
448  intvec *res;
449  int c, cp, r, t;
450 
451  t = kdim;
452  c = 1;
453  for (r=1;r<=dimtr;r++)
454  {
455  while (IMATELEM(*imat,r,c)==0) c++;
456  (*perm)[r] = c;
457  c++;
458  }
459  c = d;
460  for (r=dimtr;r>0;r--)
461  {
462  cp = (*perm)[r];
463  if (cp!=c)
464  {
465  ivKernFromRow(kern, imat, perm, t, r, c);
466  t -= (c-cp);
467  if (t==0)
468  break;
469  c = cp-1;
470  }
471  else
472  c--;
473  }
474  if (kdim>1)
475  res = ivOptimizeKern(kern);
476  else
477  res = ivTranp(kern);
478  delete kern;
479  delete perm;
480  return res;
481 }
static intvec * ivOptimizeKern(intvec *)
Definition: intvec.cc:669
static void ivKernFromRow(intvec *, intvec *, intvec *, int, int, int)
Definition: intvec.cc:613
intvec * ivTranp(intvec *o)
Definition: intvec.cc:327

◆ ivSub()

intvec* ivSub ( intvec a,
intvec b 
)

Definition at line 297 of file intvec.cc.

298 {
299  intvec * iv;
300  int mn, ma, i;
301  if (a->cols() != b->cols()) return NULL;
302  mn = si_min(a->rows(),b->rows());
303  ma = si_max(a->rows(),b->rows());
304  if (a->cols() == 1)
305  {
306  iv = new intvec(ma);
307  for (i=0; i<mn; i++) (*iv)[i] = (*a)[i] - (*b)[i];
308  if (ma > mn)
309  {
310  if (ma == a->rows())
311  {
312  for(i=mn; i<ma; i++) (*iv)[i] = (*a)[i];
313  }
314  else
315  {
316  for(i=mn; i<ma; i++) (*iv)[i] = -(*b)[i];
317  }
318  }
319  return iv;
320  }
321  if (mn != ma) return NULL;
322  iv = new intvec(a);
323  for (i=0; i<mn*a->cols(); i++) { (*iv)[i] -= (*b)[i]; }
324  return iv;
325 }

◆ ivTrace()

int ivTrace ( intvec o)

Definition at line 339 of file intvec.cc.

340 {
341  int i, s = 0, m = si_min(o->rows(),o->cols()), c = o->cols();
342  for (i=0; i<m; i++)
343  {
344  s += (*o)[i*c+i];
345  }
346  return s;
347 }

◆ ivTranp()

intvec* ivTranp ( intvec o)

Definition at line 327 of file intvec.cc.

328 {
329  int i, j, r = o->rows(), c = o->cols();
330  intvec * iv= new intvec(c, r, 0);
331  for (i=0; i<r; i++)
332  {
333  for (j=0; j<c; j++)
334  (*iv)[j*r+i] = (*o)[i*c+j];
335  }
336  return iv;
337 }

◆ ivTriangIntern()

void ivTriangIntern ( intvec imat,
int &  ready,
int &  all 
)

Definition at line 404 of file intvec.cc.

405 {
406  int rpiv, colpos=0, rowpos=0;
407  int ia=ready, ie=all;
408 
409  do
410  {
411  rowpos++;
412  do
413  {
414  colpos++;
415  rpiv = ivColPivot(imat, colpos, rowpos, ia, ie);
416  } while (rpiv==0);
417  if (rpiv>ia)
418  {
419  if (rowpos!=rpiv)
420  {
421  ivSaveRow(imat, rpiv);
422  ivFreeRow(imat, rowpos, rpiv);
423  ivSetRow(imat, rowpos, colpos);
424  rpiv = rowpos;
425  }
426  ia++;
427  if (ia==imat->cols())
428  {
429  ready = ia;
430  all = ie;
431  return;
432  }
433  }
434  ivReduce(imat, rpiv, colpos, ia, ie);
435  ivZeroElim(imat, colpos, ia, ie);
436  } while (ie>ia);
437  ready = ia;
438  all = ie;
439 }
static int ivColPivot(intvec *, int, int, int, int)
Definition: intvec.cc:484
static void ivSetRow(intvec *, int, int)
Definition: intvec.cc:513
static void ivFreeRow(intvec *, int, int)
Definition: intvec.cc:522
static void ivZeroElim(intvec *, int, int, int &)
Definition: intvec.cc:564
static void ivSaveRow(intvec *, int)
Definition: intvec.cc:505
static void ivReduce(intvec *, int, int, int, int)
Definition: intvec.cc:533

◆ ivZeroElim()

static void ivZeroElim ( intvec imat,
int  colpos,
int  ready,
int &  all 
)
static

Definition at line 564 of file intvec.cc.

566 {
567  int j, i, k, l;
568 
569  k = ready;
570  for (j=ready+1;j<=all;j++)
571  {
572  for (i=imat->cols();i>colpos;i--)
573  {
574  if (IMATELEM(*imat,j,i)!=0)
575  {
576  k++;
577  if (k<j)
578  {
579  for (l=imat->cols();l>colpos;l--)
580  IMATELEM(*imat,k,l) = IMATELEM(*imat,j,l);
581  }
582  break;
583  }
584  }
585  }
586  all = k;
587 }