$SPAD/src/input eigen.input

Timothy Daly December 3, 2016

Abstract

Contents

2

1 License

— license —

--Copyright The Numerical Algorithms Group Limited 1991.

———-

— * —

)set break resume
)spool eigen.output
)set message test on
)set message auto off
)clear all

———-

2 Computation of characteristic polynomial

2.1 For matrix of integers

— * —

--S 1 of 36
m:=matrix([[1,2,1],[2,1,-2],[1,-2,4]])
--R
--R
--R        +1   2    1 +
--R        |           |
--R   (1)  |2   1   - 2|
--R        |           |
--R        +1  - 2   4 +
--R                                                        Type: Matrix(Integer)
--E 1
--S 2 of 36
characteristicPolynomial m
--R
--R
--R            3      2
--R   (2)  - %A  + 6%A  - 25
--R                                                    Type: Polynomial(Integer)
--E 2
--S 3 of 36
characteristicPolynomial(m,x)

--R
--R
--R           3     2
--R   (3)  - x  + 6x  - 25
--R                                                    Type: Polynomial(Integer)
--E 3

———-

2.2 For matrix of polynomials

— * —

--S 4 of 36
p:=matrix([[x+1,2-x*y,x**2+1],[2-x,y+2*x,x**2-2],[y**2,x-2,4-x*y]])
--R
--R
--R        +                      2      +
--R        | x + 1   - x y + 2   x  + 1  |
--R        |                             |
--R   (4)  |                      2      |
--R        |- x + 2   y + 2x     x  - 2  |
--R        |                             |
--R        |   2                         |
--R        +  y        x - 2    - x y + 4+
--R                                            Type: Matrix(Polynomial(Integer))
--E 4
--S 5 of 36
characteristicPolynomial p
--R
--R
--R   (5)
--R         3    2           3       3            2                       2
--R     (- x  - x  + 2x - 1)y  + (- x  + (%B - 1)x  + (%B - 3)x + %B - 4)y
--R   +
--R          3             2        2                  2                 4
--R     (- 2x  + (4%B - 8)x  + (- %B  - 2%B + 16)x + %B  - 5%B + 4)y - 2x
--R   +
--R              3               2       2                   3      2
--R     (%B + 5)x  + (- 4%B + 7)x  + (3%B  - 18%B + 18)x - %B  + 5%B  + 4%B - 24
--R                                                    Type: Polynomial(Integer)
--E 5
--S 6 of 36
characteristicPolynomial(p,t)
--R
--R

--R   (6)
--R         3    2           3       3           2                     2
--R     (- x  - x  + 2x - 1)y  + (- x  + (t - 1)x  + (t - 3)x + t - 4)y
--R   +
--R          3            2       2                2                4           3
--R     (- 2x  + (4t - 8)x  + (- t  - 2t + 16)x + t  - 5t + 4)y - 2x  + (t + 5)x
--R   +
--R                2      2                 3     2
--R     (- 4t + 7)x  + (3t  - 18t + 18)x - t  + 5t  + 4t - 24
--R                                                    Type: Polynomial(Integer)
--E 6

———-

2.3 For general matrix of FRAC POLY INT

— * —

--S 7 of 36
n:=matrix([[a,b,c],[d,e,f],[g,h,k]])
--R
--R
--R        +a  b  c+
--R        |       |
--R   (7)  |d  e  f|
--R        |       |
--R        +g  h  k+
--R                                            Type: Matrix(Polynomial(Integer))
--E 7
--S 8 of 36
characteristicPolynomial n
--R
--R
--R   (8)
--R                                 2
--R     ((a - %C)e - b d - %C a + %C )k + ((- a + %C)f + c d)h
--R   +
--R                                       2                2      3
--R     (b f - c e + %C c)g + (- %C a + %C )e + %C b d + %C a - %C
--R                                                    Type: Polynomial(Integer)
--E 8

———-

There are many functions for the computation of eigenvalues and eigenvectors. The function eigenvalues returns the eigenvalues :

— * —

--S 9 of 36
leig := eigenvalues m
--R
--R
--R                  2
--R   (9)  [5,%D | %D  - %D - 5]
--IType: List(Union(Fraction(Polynomial(Integer)),SuchThat(...
--E 9
--S 10 of 36
alpha:=leig.1
--R
--R
--R   (10)  5
--R                               Type: Union(Fraction(Polynomial(Integer)),...)
--E 10

———-

alpha is a rational eigenvalue; the corresponding eigenvector can be computed with the

function eigenvector:

— * —

--S 11 of 36
eigenvector(alpha,m)
--R
--R
--R          + 0 +
--R          |   |
--R          |  1|
--R   (11)  [|- -|]
--R          |  2|
--R          |   |
--R          + 1 +
--R                  Type: List(Matrix(Fraction(Polynomial(Fraction(Integer)))))
--E 11
--S 12 of 36
beta:=leig.2
--R
--R
--R                2
--R   (12)  %D | %D  - %D - 5
--R                        Type: Union(SuchThat(Symbol,Polynomial(Integer)),...)
--E 12

———-

beta is a ”symbolic” eigenvalue, i.e. it is a root (symbolically expressed) of an irreducible factor of the characteristic polynomial. In this case too we can compute the associate eigenvectors.

— * —

--S 13 of 36
eigenvector(beta,m)$EP(INT)
--R
--R
--R          +%D+
--R          |  |
--R   (13)  [|2 |]
--R          |  |
--R          +1 +
--R                            Type: List(Matrix(Fraction(Polynomial(Integer))))
--E 13
-- eigenvector(beta,m)  not accepted by the interpreter

———-

The function eigenvectors can compute Eigenvalues and eigenvectors simultaneously with the function

— * —

--S 14 of 36
eigenvectors m
--R
--R
--R   (14)
--R                                   + 0 +
--R                                   |   |
--R                                   |  1|
--R   [[eigval= 5,eigmult= 1,eigvec= [|- -|]],
--R                                   |  2|
--R                                   |   |
--R                                   + 1 +
--R                                                     +%E+
--R                     2                               |  |
--R    [eigval= (%E | %E  - %E - 5),eigmult= 1,eigvec= [|2 |]]]
--R                                                     |  |
--R                                                     +1 +
--IType: List(Record(eigval: Union(Fraction(Polynomial(Integer)),...
--E 14
--S 15 of 36
q:=matrix [[x**2-y**2,(x-y)*(2*x+3*y)],[x+y,2*x+3*y]]
--R
--R
--R         +   2    2      2           2+
--R   (15)  |- y  + x   - 3y  + x y + 2x |
--R         |                            |

--R         +  y + x         3y + 2x     +
--R                                            Type: Matrix(Polynomial(Integer))
--E 15
--S 16 of 36
eigenvectors(q)
--R
--R
--R   (16)
--R                2         2                          +- y + x+
--R   [[eigval= - y  + 3y + x  + 2x,eigmult= 1,eigvec= [|       |]],
--R                                                     +   1   +
--R                                   +- 3y - 2x+
--R                                   |---------|
--R    [eigval= 0,eigmult= 1,eigvec= [|  y + x  |]]]
--R                                   |         |
--R                                   +    1    +
--IType: List(Record(eigval: Union(Fraction(Polynomial(Integer)),...
--E 16
--S 17 of 36
p:=matrix([[76,-18,58,-10],[-4,78,2,-2],[-6,15,45,3],[22,-75,7,41]])
--R
--R
--R         +76   - 18  58  - 10+
--R         |                   |
--R         |- 4   78   2   - 2 |
--R   (17)  |                   |
--R         |- 6   15   45   3  |
--R         |                   |
--R         +22   - 75  7    41 +
--R                                                        Type: Matrix(Integer)
--E 17
--S 18 of 36
ll := eigenvectors p
--R
--R
--R   (18)
--R                                    +10 +
--R                                    |-- |
--R                                    | 7 |
--R                                    |   |
--R                                    | 2 |
--R                                    | - |
--R   [[eigval= 48,eigmult= 2,eigvec= [| 7 |]],
--R                                    |   |
--R                                    |  3|
--R                                    |- -|
--R                                    |  7|

--R                                    |   |
--R                                    + 1 +
--R                                    +- 2+
--R                                    |   |
--R                                    |- 1|
--R    [eigval= 72,eigmult= 2,eigvec= [|   |]]]
--R                                    | 0 |
--R                                    |   |
--R                                    + 1 +
--IType: List(Record(eigval: Union(Fraction(Polynomial(Integer)),...
--E 18

———-

In this case the algebraic multiplicity (the field eigmult) is different from the geometric multiplicity (the length of the field eigvec).

— * —

--S 19 of 36
generalizedEigenvectors p
--R
--R
--R   (19)
--R                            +  10+
--R                            |- --|
--R                            |   3| +0+                           +- 12+ +- 2+
--R                            |    | | |                           |    | |   |
--R                            |  2 | |0|                           |- 3 | |- 1|
--R   [[eigval= 48,geneigvec= [|- - |,| |]],[eigval= 72,geneigvec= [|    |,|   |]]]
--R                            |  3 | |0|                           | 1  | | 0 |
--R                            |    | | |                           |    | |   |
--R                            | 1  | +1+                           + 0  + + 1 +
--R                            |    |
--R                            + 0  +
--IType: List(Record(eigval: Union(Fraction(Polynomial(Integer)),...
--E 19
--S 20 of 36
generalizedEigenvector(ll.1,p)$EP(INT)
--R
--R
--R          +  10+
--R          |- --|
--R          |   3| +0+
--R          |    | | |
--R          |  2 | |0|
--R   (20)  [|- - |,| |]
--R          |  3 | |0|
--R          |    | | |

--R          | 1  | +1+
--R          |    |
--R          + 0  +
--R                            Type: List(Matrix(Fraction(Polynomial(Integer))))
--E 20
-- generalizedEigenvector(ll.1,p) the interpreter can not handle this

———-

These functions return respectively the complete set of generalized eigenvectors or the generalized eigenvectors associated to a particular eigenvalue alpha, i.e. a basis of the nullSpace((palpha*I)**k) where k is the algebraic multiplicity of alpha.

In the case of symbolic eigenvalues it is possible to convert the symbolic eigenvalue and the corresponding eigenvectors in a more explicit form.

— * —

--S 21 of 36
m
--R
--R
--R         +1   2    1 +
--R         |           |
--R   (21)  |2   1   - 2|
--R         |           |
--R         +1  - 2   4 +
--R                                                        Type: Matrix(Integer)
--E 21
--S 22 of 36
mm:=matrix([[30,4,24],[-17,8,-7],[-31,-54,-5]])
--R
--R
--R         + 30    4    24 +
--R         |               |
--R   (22)  |- 17   8    - 7|
--R         |               |
--R         +- 31  - 54  - 5+
--R                                                        Type: Matrix(Integer)
--E 22

———-

The function radicalEigenvalues expresses, when possible, the eigenvalues in terms of radicals.

— * —

--S 23 of 36
le1:=radicalEigenvalues m
--R

--R
--R             +--+      +--+
--R          - \|21  + 1 \|21  + 1
--R   (23)  [-----------,---------,5]
--R               2          2
--R                                              Type: List(Expression(Integer))
--E 23
--S 24 of 36
le2:=radicalEigenvalues mm
--R
--R
--R             +---+      +---+
--R   (24)  [22\|- 1 ,- 22\|- 1 ,33]
--R                                              Type: List(Expression(Integer))
--E 24

———-

The function radicalEigenvector computes the eigenvectors associated to a given eigenvalue, expressed in terms of radicals

— * —

--S 25 of 36
radicalEigenvector(le1.2, m)
--R
--R
--R          +    10   +
--R          |---------|
--R          | +--+    |
--R   (25)  [|\|21  - 1|]
--R          |         |
--R          |    2    |
--R          |         |
--R          +    1    +
--R                                      Type: List(Matrix(Expression(Integer)))
--E 25
--S 26 of 36
radicalEigenvector(le2.2,mm)
--R
--R
--R          +       +---+       +
--R          |- 1449\|- 1  + 1720|
--R          |-------------------|
--R          |      +---+        |
--R          |  328\|- 1  - 3343 |
--R          |                   |
--R   (26)  [|      +---+        |]
--R          |    7\|- 1  - 9    |

--R          |   ------------    |
--R          |      +---+        |
--R          |   38\|- 1  - 8    |
--R          |                   |
--R          +         1         +
--R                                      Type: List(Matrix(Expression(Integer)))
--E 26

———-

The function radicalEigenvectors computes simultaneously all the eigenvalues and the associated eigenvectors and expresses, when it is possible, in terms of radical.

— * —

--S 27 of 36
radicalEigenvectors m
--R
--R
--R   (27)
--R                                            + +--+    +
--R              +--+                          |\|21  + 1|
--R             \|21  + 1                      |---------|
--R   [[radval= ---------,radmult= 1,radvect= [|    2    |]],
--R                 2                          |         |
--R                                            |    2    |
--R                                            |         |
--R                                            +    1    +
--R                                              +   +--+    +
--R                +--+                          |- \|21  + 1|
--R             - \|21  + 1                      |-----------|
--R    [radval= -----------,radmult= 1,radvect= [|     2     |]],
--R                  2                           |           |
--R                                              |     2     |
--R                                              |           |
--R                                              +     1     +
--R                                    + 0 +
--R                                    |   |
--R                                    |  1|
--R    [radval= 5,radmult= 1,radvect= [|- -|]]]
--R                                    |  2|
--R                                    |   |
--R                                    + 1 +
--IType: List(Record(radval: Expression(Integer),radmult: Integer,...
--E 27
--S 28 of 36
radicalEigenvectors mm
--R
--R
--R   (28)

--R                                             +   +---+     +
--R                                             |11\|- 1  - 16|
--R                                             |-------------|
--R                                             |      29     |
--R                  +---+                      |             |
--R   [[radval= - 22\|- 1 ,radmult= 1,radvect= [|   +---+     |]],
--R                                             |11\|- 1  + 13|
--R                                             |-------------|
--R                                             |      58     |
--R                                             |             |
--R                                             +      1      +
--R                                           +     +---+     +
--R                                           |- 11\|- 1  - 16|
--R                                           |---------------|
--R                                           |       29      |
--R                +---+                      |               |
--R    [radval= 22\|- 1 ,radmult= 1,radvect= [|     +---+     |]],
--R                                           |- 11\|- 1  + 13|
--R                                           |---------------|
--R                                           |       58      |
--R                                           |               |
--R                                           +       1       +
--R                                     + 4 +
--R                                     |   |
--R    [radval= 33,radmult= 1,radvect= [|- 3|]]]
--R                                     |   |
--R                                     + 1 +
--IType: List(Record(radval: Expression(Integer),radmult: Integer,...
--E 28

———-

There exist analogous functions for the computation of real and complex eigenvalues and eigenvectors.

In order to compute respectively the eigenvalues or the eigenvalues and the associared eigevectors of a matrix m and express them as rational numbers (Gaussian rational) up to precision 1/1000000 use

— * —

--S 29 of 36
realEigenvalues(m,1/1000000)
--R
--R
--R            3756603   5853755
--R   (29)  [- -------,5,-------]
--R            2097152   2097152
--R                                                Type: List(Fraction(Integer))
--E 29
--S 30 of 36

complexEigenvalues(mm,1/1000000)
--R
--R
--R   (30)  [- 22%i,22%i,33]
--R                                       Type: List(Complex(Fraction(Integer)))
--E 30
--S 31 of 36
realEigenvectors(m,1/1000000)
--R
--R
--R   (31)
--R                                    + 0 +
--R                                    |   |
--R                                    |  1|
--R   [[outval= 5,outmult= 1,outvect= [|- -|]],
--R                                    |  2|
--R                                    |   |
--R                                    + 1 +
--R                                          +5853755+
--R                                          |-------|
--R             5853755                      |2097152|
--R    [outval= -------,outmult= 1,outvect= [|       |]],
--R             2097152                      |   2   |
--R                                          |       |
--R                                          +   1   +
--R                                            +  3756603+
--R                                            |- -------|
--R               3756603                      |  2097152|
--R    [outval= - -------,outmult= 1,outvect= [|         |]]]
--R               2097152                      |    2    |
--R                                            |         |
--R                                            +    1    +
--IType: List(Record(outval: Fraction(Integer),outmult: Integer,...
--E 31
--S 32 of 36
complexEigenvectors(mm,1/1000000)
--R
--R
--R   (32)
--R                                     + 4 +
--R                                     |   |
--R   [[outval= 33,outmult= 1,outvect= [|- 3|]],
--R                                     |   |
--R                                     + 1 +
--R                                       +  16   11   +
--R                                       |- -- - -- %i|
--R                                       |  29   29   |
--R                                       |            |

--R    [outval= 22%i,outmult= 1,outvect= [| 13   11    |]],
--R                                       | -- - -- %i |
--R                                       | 58   58    |
--R                                       |            |
--R                                       +     1      +
--R                                         +  16   11   +
--R                                         |- -- + -- %i|
--R                                         |  29   29   |
--R                                         |            |
--R    [outval= - 22%i,outmult= 1,outvect= [| 13   11    |]]]
--R                                         | -- + -- %i |
--R                                         | 58   58    |
--R                                         |            |
--R                                         +     1      +
--IType: List(Record(outval: Complex(Fraction(Integer)),...
--E 32

———-

To have the eigenvalues expressed as real float (complex float) use realEigenvalues — * —

--S 33 of 36
realEigenvalues(m,.000001)
--R
--R
--R   (33)  [- 1.7912878990 173339844,5.0,2.7912878990 173339844]
--R                                                            Type: List(Float)
--E 33
--S 34 of 36
realEigenvectors(m,.000001)
--R
--R
--R   (34)
--R                                      + 0.0 +
--R                                      |     |
--R   [[outval= 5.0,outmult= 1,outvect= [|- 0.5|]],
--R                                      |     |
--R                                      + 1.0 +
--R
--R     [outval= 2.7912878990 173339844, outmult= 1,
--R                +2.7912878990 173339844+
--R                |                      |
--R      outvect= [|         2.0          |]]
--R                |                      |
--R                +         1.0          +
--R     ,
--R
--R     [outval= - 1.7912878990 173339844, outmult= 1,

--R                +- 1.7912878990 173339844+
--R                |                        |
--R      outvect= [|          2.0           |]]
--R                |                        |
--R                +          1.0           +
--R     ]
--RType: List(Record(outval: Float,outmult: Integer,outvect: List(Matrix(Float))))
--E 34
--S 35 of 36
complexEigenvalues(mm,.000001)
--R
--R
--R   (35)  [- 22.0 %i,22.0 %i,33.0]
--R                                                   Type: List(Complex(Float))
--E 35
--S 36 of 36
complexEigenvectors(mm,.000001)
--R
--R
--R   (36)
--R                                       + 4.0 +
--R                                       |     |
--R   [[outval= 33.0,outmult= 1,outvect= [|- 3.0|]],
--R                                       |     |
--R                                       + 1.0 +
--R
--R     [outval= 22.0 %i, outmult= 1,
--R                +- 0.5517241379 3103448276 - 0.3793103448 275862069 %i+
--R                |                                                     |
--R      outvect= [|0.2241379310 3448275862 - 0.1896551724 1379310345 %i |]]
--R                |                                                     |
--R                +                         1.0                         +
--R     ,
--R
--R     [outval= - 22.0 %i, outmult= 1,
--R                +- 0.5517241379 3103448276 + 0.3793103448 275862069 %i+
--R                |                                                     |
--R      outvect= [|0.2241379310 3448275862 + 0.1896551724 1379310345 %i |]]
--R                |                                                     |
--R                +                         1.0                         +
--R     ]
--IType: List(Record(outval: Complex(Float),outmult: Integer,...
--E 36
)spool
)lisp (bye)

———-

References^{[1] nothing}

17