Download this file

EigenLib3.def    221 lines (198 with data), 14.2 kB

  1
  2
  3
  4
  5
  6
  7
  8
  9
 10
 11
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
DEFINITION MODULE EigenLib3;
(*------------------------------------------------------------------------*)
(* Routinen zu Berechnung der Eigensysteme komplexer genereller Matrizen. *)
(* Routines to calculate the eigensystems of complex general matrices. *)
(*------------------------------------------------------------------------*)
(* Implementation : Michael Riedl *)
(* Licence : GNU Lesser General Public License (LGPL) *)
(*------------------------------------------------------------------------*)
(* $Id: EigenLib3.def,v 1.4 2017/10/22 17:51:58 mriedl Exp mriedl $ *)
FROM Deklera IMPORT CARDVEKTOR,CVEKTOR,CMATRIX;
PROCEDURE CSortEwEv(VAR EW : CVEKTOR; (* Eigenwerte *)
VAR EV : CMATRIX; (* Eigenvektoren *)
dim : CARDINAL;
ord : INTEGER);
(*----------------------------------------------------------*)
(* Sortiert Eigenwerte und Eigenvektoren nach Gr"o\3e der *)
(* Absolutwerte der Eigenwerte. *)
(* Wenn ord = 1 in absteigender Reihenfolge, sonst in *)
(* ansteigender Reihenfolge. *)
(*----------------------------------------------------------*)
PROCEDURE CBalance(VAR A : ARRAY OF ARRAY OF LONGCOMPLEX;
dim : CARDINAL;
VAR low,high : CARDINAL;
VAR Scale : ARRAY OF LONGREAL;
VAR iFehl : INTEGER);
(*----------------------------------------------------------------*)
(* Balanciert eine quadratische complexe Matrix als Vorbereitung *)
(* zur Eigensystemberechnung) *)
(* Werden die Eigenvektoren von A ben"otigt, ist nach deren *)
(* Berechnung die Prozedur BalBak anzuwenden ! *)
(* *)
(* Reduce the norm of a [1:n,1:n] by exact diagonal similarity *)
(* transformations stored in Scale[1:n] *)
(* *)
(* Lit.: Parlett,B.N.; Reinsch,C.; Numer. Math 13, 293 (1969) *)
(* Handbook - Algolprozedur cbalance. *)
(*----------------------------------------------------------------*)
PROCEDURE CBalBak( dim : CARDINAL; (* Ordnung der Matrix Z *)
low,high : CARDINAL;
VAR Scale : ARRAY OF LONGREAL;
m : CARDINAL;
VAR Z : ARRAY OF ARRAY OF LONGCOMPLEX);
(*----------------------------------------------------------------*)
(* Formt die Eigenvektoren einer complxen, unsymmetrischen Matrix *)
(* durch R"ucktransformation der Eigenvektoren einer durch die *)
(* Prozedur CBalance 'balancierten' Matrix. *)
(* *)
(* ==> Scale : enth"alt Informationen "uber die Permutationen *)
(* und Skalierungsfaktoren, die von der Prozedur *)
(* Balance benutzt wurden. *)
(* m : Anzahl der Vektoren in Z, die R"ucktrans- *)
(* formiert werden sollen. *)
(* <==> Z : Eigenvektoren, die R"ucktransformiert werden *)
(* sollen (in den ersten m Zeilen). *)
(* *)
(* Backward transformation of a set of right-hand eigenvectors *)
(* of a balanced matrix into the eigenvectors of the original *)
(* matrix from which the balanced matrix was derived by a call of *)
(* procedure balance *)
(* *)
(* Lit.: Parlett; Reinsch; Num. Math. 13, 293-304 (1969) *)
(* Handbook - Algolprozedur cbalbak. *)
(*----------------------------------------------------------------*)
PROCEDURE ComHess(VAR A : CMATRIX;
dim : CARDINAL;
k,l : CARDINAL;
VAR Index : CARDVEKTOR);
(*----------------------------------------------------------------*)
(* Reduziert die in A "ubergebene, unsymmetrische komplexe Matrix *)
(* auf obere Hessenbergform. *)
(* Unterhalb der Haupdiagonalen und in dem Feld Index werden die *)
(* Informationen "uber die Transformationsschritte gespeichert. *)
(* Die Parameter k und l definieren die Submatrix, die trans- *)
(* formiert werden soll (A[k,k] bis A[l,l]). *)
(* Siehe auch : ComHessBak. *)
(* Lit. : Martin, R.S. et. al.; Numer. Math. 12, 349 (1968) *)
(* Handbook - Algolprozedur comhes. *)
(*----------------------------------------------------------------*)
PROCEDURE ComHessBak(VAR EV : CMATRIX;
VAR A : CMATRIX;
Index : CARDVEKTOR;
dim : CARDINAL;
EvAnz : CARDINAL;
k,l : CARDINAL);
(*----------------------------------------------------------------*)
(* Berechent aus den Eigenvektoren EV der in A gespeicherten *)
(* oberen Hessenbergmatrix die Eigenvektoren der urspr"unglichen, *)
(* mit der Prozedur ComHess reduzierten Matrix. *)
(* Dazu mu\3 A und das Feld Index in der Form, wie es in ComHess *)
(* berechnet wurde, unver"andert "ubergeben werden. *)
(* EV : Eingang : Eigenvektoren der Hessenbergmatrix. *)
(* Ausgang : Eigenvektoren der urspr"unglichen Matrix. *)
(* Index : Transformationsinformationen, in ComHess berechent. *)
(* dim : Aktuelle Dimension von A. *)
(* EvAnz : Alle Eigenvektoren von 1..EvAnz werden transformiert *)
(* k,l : Definiert Submatrix, die berechent wurde. *)
(* *)
(* Lit. : Martin, R.S. et. al.; Numer. Math. 12, 349 (1968) *)
(* Handbook - Algolprozedur combak. *)
(*----------------------------------------------------------------*)
PROCEDURE CHessLR(VAR H : CMATRIX; (* komplexe Hessenbergmatrix *)
VAR EW : CVEKTOR; (* komplexe Eigenwerte von H *)
dim : CARDINAL); (* Dimension von H *)
(*----------------------------------------------------------------*)
(* Berechnet die Eigenwerte der komplexen oberen Hessenberg- *)
(* matrix H. Die Matrix H wird dabei zerstoert. *)
(* Lit.: Matrin,R.S.; Wilkinson,J.H.; Numer. Math. 12, 369 (1968) *)
(* Handbook - Algolprozedure comlr. *)
(*----------------------------------------------------------------*)
PROCEDURE CHessInvIter(VAR EV : CMATRIX; (* Eigenvektoren *)
VAR EW : CVEKTOR; (* Eigenwerte von A *)
VAR H : CMATRIX; (* Hessenbergmatrix *)
dim : CARDINAL; (* Dimension von H *)
VAR IErr : ARRAY OF CARDINAL;
VAR IFehl : INTEGER);
(*----------------------------------------------------------------*)
(* Berechnet die Eigenvektoren einer komplexen oberen Hessen- *)
(* bergmatrix H. *)
(* *)
(* Given a complex upper Hessenberg matrix in the arrays *)
(* H[1:n,1:n] and its eigenvalues in EW[1:n] this procedure finds *)
(* the eigen-vectors and stores them in arrays EV[1:n,1:n]. Any *)
(* vector which has not been accepted is set to zero. *)
(* *)
(* Input *)
(* *)
(* c : an nx1 boolean array indicating the eigenvalues for *)
(* which the eigenvectors are required (NOT USED) *)
(* r : as used below is the number of components which are *)
(* true (NOT USED) *)
(* H : nxn complex arrays giving the matrix *)
(* EW : nx1 complex arrays giving the real & imaginary parts *)
(* of the eigenvalues of A. These must be as given by *)
(* comlr since eigenvalues must be ordered to have the *)
(* correct affiliations *)
(* *)
(* Output *)
(* *)
(* EW : the real parts of the selected eigenvalues, *)
(* perturbed where necessary in the case of close *)
(* eigenvalues. *)
(* EV : two nxr arrays giving the real and imaginary parts *)
(* of the r eigenvectors; *)
(* IErr : IErr[i] = k : Nummer des nicht berechneten Eigen- *)
(* vektors. ( IErr[1..IFehl] ) *)
(* IFehl : Anzahl der fehlgeschlagenen Vektoriterationen. *)
(* *)
(* *)
(* Lit.: Peters,G.; Wilkinson, J. H.: Eigenvectors of real and *)
(* complex matrices by LR and QR triangularizations. *)
(* Numer. Math. 16, 181-204 (1970). *)
(* Handbook - Algolprozedur cx invit. Cf. II/15 *)
(*----------------------------------------------------------------*)
PROCEDURE ComEigEvEw(VAR A : CMATRIX;
dim : CARDINAL;
VAR EW : CVEKTOR;
VAR EV : CMATRIX;
Norm : CARDINAL;
ord : INTEGER;
VAR iFehl : INTEGER);
(*----------------------------------------------------------------*)
(* Berechnet die Eigenvektoren einer komplexen generellen Matrix *)
(* A mithilfe der Prozeduren ComHess,ComHessBak,CHessLR und *)
(* CHessInvIter *)
(* *)
(* ==> A : Komplexe Eingabematrix *)
(* ==> dim : Dimension der Matrix A *)
(* <== EW : Eigenwerte *)
(* <== EV : Rechtsvektoren *)
(* ==> Norm : Wenn 0 keine Normierung, wenn 1 Maximumnorm, *)
(* wenn 2 Euclidsche Norm *)
(* ==> ord : siehe Prozedur CSortEwEv *)
(* <==> iFehl : siehe Prozedur CHessInvIter *)
(*----------------------------------------------------------------*)
PROCEDURE ComEig(VAR A : ARRAY OF ARRAY OF LONGCOMPLEX;
VAR REV : ARRAY OF ARRAY OF LONGCOMPLEX; (* Rechts-EV *)
VAR LEV : ARRAY OF ARRAY OF LONGCOMPLEX; (* Links-EV *)
VAR EW : ARRAY OF LONGCOMPLEX; (* Eigenwerte *)
dim : CARDINAL;
genau : LONGREAL;
ord : INTEGER);
(*----------------------------------------------------------------*)
(* Berechnet die Eigenvektoren einer komplexen generellen *)
(* matrix A nach der generalisierten Jacobi-Methode *)
(* *)
(* ==> A : Komplexe Eingabematrix *)
(* <== REV : Rechtsvektoren *)
(* <== LEV : Linksvektoren *)
(* <== EW : Eigenwerte *)
(* ==> dim : Dimension der Matrix A *)
(* ==> genau : Vorgabewert ab wann der Absultwert eines nicht- *)
(* diagonalelements als Null angesehen werden soll *)
(* ==> ord : siehe Prozedur CSortEwEv *)
(* *)
(* Lit. : P.J.Eberlein, Numer. Math. 14, 232 (1970) *)
(* Handbook - Algolprozedur comeig. *)
(*----------------------------------------------------------------*)
END EigenLib3.