eigenval.cc
Go to the documentation of this file.
1 /*****************************************
2 * Computer Algebra System SINGULAR *
3 *****************************************/
4 /*
5 * ABSTRACT: eigenvalues of constant square matrices
6 */
7 
8 
9 
10 
11 #include <kernel/mod2.h>
12 
13 #ifdef HAVE_EIGENVAL
14 
15 #include <kernel/structs.h>
16 #include <misc/intvec.h>
17 #include <coeffs/numbers.h>
18 #include <kernel/polys.h>
19 #include <kernel/ideals.h>
20 #include <polys/matpol.h>
21 #include <polys/clapsing.h>
23 
24 
25 matrix evSwap(matrix M,int i,int j)
26 {
27  if(i==j)
28  return(M);
29 
30  for(int k=1;k<=MATROWS(M);k++)
31  {
32  poly p=MATELEM(M,i,k);
33  MATELEM(M,i,k)=MATELEM(M,j,k);
34  MATELEM(M,j,k)=p;
35  }
36 
37  for(int k=1;k<=MATCOLS(M);k++)
38  {
39  poly p=MATELEM(M,k,i);
40  MATELEM(M,k,i)=MATELEM(M,k,j);
41  MATELEM(M,k,j)=p;
42  }
43 
44  return(M);
45 }
46 
47 matrix evRowElim(matrix M,int i,int j,int k)
48 {
49  if(MATELEM(M,i,k)==NULL||MATELEM(M,j,k)==NULL)
50  return(M);
51  poly p1=pp_Jet(MATELEM(M,i,k),0,currRing);
52  poly p2=pp_Jet(MATELEM(M,j,k),0,currRing);
53  if ((p1==NULL)||(p2==NULL)) return (M);
54 
55  poly p=pNSet(nDiv(pGetCoeff(p1),pGetCoeff(p2)));
56  pNormalize(p);
57 
58  for(int l=1;l<=MATCOLS(M);l++)
59  {
60  MATELEM(M,i,l)=pSub(MATELEM(M,i,l),ppMult_qq(p,MATELEM(M,j,l)));
61  pNormalize(MATELEM(M,i,l));
62  }
63  for(int l=1;l<=MATROWS(M);l++)
64  {
65  MATELEM(M,l,j)=pAdd(MATELEM(M,l,j),ppMult_qq(p,MATELEM(M,l,i)));
66  pNormalize(MATELEM(M,l,j));
67  }
68 
69  pDelete(&p);
70  pDelete(&p1);
71  pDelete(&p2);
72 
73  return(M);
74 }
75 
76 matrix evColElim(matrix M,int i,int j,int k)
77 {
78  if(MATELEM(M,k,i)==0||MATELEM(M,k,j)==0)
79  return(M);
80 
81  poly p=pNSet(nDiv(pGetCoeff(MATELEM(M,k,i)),pGetCoeff(MATELEM(M,k,j))));
82  pNormalize(p);
83 
84  for(int l=1;l<=MATROWS(M);l++)
85  {
86  MATELEM(M,l,i)=pSub(MATELEM(M,l,i),ppMult_qq(p,MATELEM(M,l,j)));
87  pNormalize(MATELEM(M,l,i));
88  }
89  for(int l=1;l<=MATCOLS(M);l++)
90  {
91  MATELEM(M,j,l)=pAdd(MATELEM(M,j,l),ppMult_qq(p,MATELEM(M,i,l)));
92  pNormalize(MATELEM(M,j,l));
93  }
94 
95  pDelete(&p);
96 
97  return(M);
98 }
99 
101 {
102  int n=MATROWS(M);
103  if(n!=MATCOLS(M))
104  return(M);
105 
106  for(int k=1,j=2;k<n-1;k++,j=k+1)
107  {
108  while((j<=n)
109  &&((MATELEM(M,j,k)==NULL)
110  || (p_Totaldegree(MATELEM(M,j,k),currRing)!=0)))
111  j++;
112 
113  if(j<=n)
114  {
115  M=evSwap(M,j,k+1);
116 
117  for(int i=j+1;i<=n;i++)
118  M=evRowElim(M,i,k+1,k);
119  }
120  }
121 
122  return(M);
123 }
124 
125 #endif /* HAVE_EIGENVAL */
#define ppMult_qq(p, q)
Definition: polys.h:191
#define pAdd(p, q)
Definition: polys.h:186
matrix evRowElim(matrix M, int i, int j, int k)
Definition: eigenval.cc:47
#define pNSet(n)
Definition: polys.h:296
Compatiblity layer for legacy polynomial operations (over currRing)
return P p
Definition: myNF.cc:203
static long p_Totaldegree(poly p, const ring r)
Definition: p_polys.h:1430
int k
Definition: cfEzgcd.cc:93
static number & pGetCoeff(poly p)
return an alias to the leading coefficient of p assumes that p != NULL NOTE: not copy ...
Definition: monomials.h:51
matrix evSwap(matrix M, int i, int j)
Definition: eigenval.cc:25
#define M
Definition: sirandom.c:24
ring currRing
Widely used global variable which specifies the current polynomial ring for Singular interpreter and ...
Definition: polys.cc:10
#define pSub(a, b)
Definition: polys.h:270
int j
Definition: myNF.cc:70
pNormalize(P.p)
int i
Definition: cfEzgcd.cc:123
matrix evHessenberg(matrix M)
Definition: eigenval.cc:100
#define nDiv(a, b)
Definition: numbers.h:32
#define MATCOLS(i)
Definition: matpol.h:28
#define NULL
Definition: omList.c:10
#define pDelete(p_ptr)
Definition: polys.h:169
END_NAMESPACE const void * p2
Definition: syzextra.cc:202
poly pp_Jet(poly p, int m, const ring R)
Definition: p_polys.cc:4167
#define MATROWS(i)
Definition: matpol.h:27
polyrec * poly
Definition: hilb.h:10
matrix evColElim(matrix M, int i, int j, int k)
Definition: eigenval.cc:76
int l
Definition: cfEzgcd.cc:94
#define MATELEM(mat, i, j)
Definition: matpol.h:29