[35aab3] | 1 | /***************************************** |
---|
| 2 | * Computer Algebra System SINGULAR * |
---|
| 3 | *****************************************/ |
---|
[341696] | 4 | /* $Id$ */ |
---|
[35aab3] | 5 | /* |
---|
| 6 | * ABSTRACT: eigenvalues of constant square matrices |
---|
| 7 | */ |
---|
| 8 | |
---|
[599326] | 9 | #include <kernel/mod2.h> |
---|
[35aab3] | 10 | |
---|
| 11 | #ifdef HAVE_EIGENVAL |
---|
| 12 | |
---|
[599326] | 13 | #include <kernel/febase.h> |
---|
| 14 | #include <kernel/structs.h> |
---|
[76b4bd] | 15 | //#include "ipid.h" |
---|
[599326] | 16 | #include <kernel/intvec.h> |
---|
| 17 | #include <kernel/numbers.h> |
---|
| 18 | #include <kernel/polys.h> |
---|
| 19 | #include <kernel/ideals.h> |
---|
[76b4bd] | 20 | //#include "lists.h" |
---|
[599326] | 21 | #include <kernel/matpol.h> |
---|
| 22 | #include <kernel/clapsing.h> |
---|
| 23 | #include <kernel/eigenval.h> |
---|
[35aab3] | 24 | |
---|
| 25 | |
---|
| 26 | matrix evSwap(matrix M,int i,int j) |
---|
| 27 | { |
---|
| 28 | if(i==j) |
---|
| 29 | return(M); |
---|
| 30 | |
---|
| 31 | for(int k=1;k<=MATROWS(M);k++) |
---|
| 32 | { |
---|
| 33 | poly p=MATELEM(M,i,k); |
---|
| 34 | MATELEM(M,i,k)=MATELEM(M,j,k); |
---|
| 35 | MATELEM(M,j,k)=p; |
---|
| 36 | } |
---|
| 37 | |
---|
| 38 | for(int k=1;k<=MATCOLS(M);k++) |
---|
| 39 | { |
---|
| 40 | poly p=MATELEM(M,k,i); |
---|
| 41 | MATELEM(M,k,i)=MATELEM(M,k,j); |
---|
| 42 | MATELEM(M,k,j)=p; |
---|
| 43 | } |
---|
| 44 | |
---|
| 45 | return(M); |
---|
| 46 | } |
---|
| 47 | |
---|
| 48 | matrix evRowElim(matrix M,int i,int j,int k) |
---|
| 49 | { |
---|
| 50 | if(MATELEM(M,i,k)==NULL||MATELEM(M,j,k)==NULL) |
---|
| 51 | return(M); |
---|
| 52 | |
---|
| 53 | poly p=pNSet(nDiv(pGetCoeff(MATELEM(M,i,k)),pGetCoeff(MATELEM(M,j,k)))); |
---|
| 54 | pNormalize(p); |
---|
| 55 | |
---|
| 56 | for(int l=1;l<=MATCOLS(M);l++) |
---|
| 57 | { |
---|
| 58 | MATELEM(M,i,l)=pSub(MATELEM(M,i,l),ppMult_qq(p,MATELEM(M,j,l))); |
---|
| 59 | pNormalize(MATELEM(M,i,l)); |
---|
| 60 | } |
---|
| 61 | for(int l=1;l<=MATROWS(M);l++) |
---|
| 62 | { |
---|
| 63 | MATELEM(M,l,j)=pAdd(MATELEM(M,l,j),ppMult_qq(p,MATELEM(M,l,i))); |
---|
| 64 | pNormalize(MATELEM(M,l,j)); |
---|
| 65 | } |
---|
| 66 | |
---|
| 67 | pDelete(&p); |
---|
| 68 | |
---|
| 69 | return(M); |
---|
| 70 | } |
---|
| 71 | |
---|
| 72 | matrix evColElim(matrix M,int i,int j,int k) |
---|
| 73 | { |
---|
| 74 | if(MATELEM(M,k,i)==0||MATELEM(M,k,j)==0) |
---|
| 75 | return(M); |
---|
| 76 | |
---|
| 77 | poly p=pNSet(nDiv(pGetCoeff(MATELEM(M,k,i)),pGetCoeff(MATELEM(M,k,j)))); |
---|
| 78 | pNormalize(p); |
---|
| 79 | |
---|
| 80 | for(int l=1;l<=MATROWS(M);l++) |
---|
| 81 | { |
---|
| 82 | MATELEM(M,l,i)=pSub(MATELEM(M,l,i),ppMult_qq(p,MATELEM(M,l,j))); |
---|
| 83 | pNormalize(MATELEM(M,l,i)); |
---|
| 84 | } |
---|
| 85 | for(int l=1;l<=MATCOLS(M);l++) |
---|
| 86 | { |
---|
| 87 | MATELEM(M,j,l)=pAdd(MATELEM(M,j,l),ppMult_qq(p,MATELEM(M,i,l))); |
---|
| 88 | pNormalize(MATELEM(M,j,l)); |
---|
| 89 | } |
---|
| 90 | |
---|
| 91 | pDelete(&p); |
---|
| 92 | |
---|
| 93 | return(M); |
---|
| 94 | } |
---|
| 95 | |
---|
| 96 | matrix evHessenberg(matrix M) |
---|
| 97 | { |
---|
| 98 | int n=MATROWS(M); |
---|
| 99 | if(n!=MATCOLS(M)) |
---|
| 100 | return(M); |
---|
| 101 | |
---|
| 102 | for(int k=1,j=2;k<n-1;k++,j=k+1) |
---|
| 103 | { |
---|
| 104 | while(j<=n&&MATELEM(M,j,k)==0) |
---|
| 105 | j++; |
---|
| 106 | |
---|
| 107 | if(j<=n) |
---|
| 108 | { |
---|
| 109 | M=evSwap(M,j,k+1); |
---|
| 110 | |
---|
| 111 | for(int i=j+1;i<=n;i++) |
---|
| 112 | M=evRowElim(M,i,k+1,k); |
---|
| 113 | } |
---|
| 114 | } |
---|
| 115 | |
---|
| 116 | return(M); |
---|
| 117 | } |
---|
| 118 | |
---|
| 119 | #endif /* HAVE_EIGENVAL */ |
---|