-
Notifications
You must be signed in to change notification settings - Fork 0
/
eignv.c
43 lines (37 loc) · 914 Bytes
/
eignv.c
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
#include <math.h>
#include "constant.h"
#include "functions.h"
void eigv (a, q, size)
double **a, **q;
int size;
{
int is, il, first=0, last, shortest, longest = size;
ident (q, size, size);
m_hes (a, q, size);
while (longest > 1) {
last=size-1, shortest=size-1, longest = 0;
shortest = size-1;
il = 0; is = 0;
while (is < size-1) {
while (is<size-1 &&
(fabs(a[is+1][is])>EPS*fabs(a[is][is]) ||
fabs(a[is+1][is])>EPS*fabs(a[is+1][is+1])))
is++; /* "is" is next zero sub-diag element */
addflops (1);
if (is<size-1)
a[is+1][is] = 0.0; /* Prevent further pivots */
if ((is-il > 1) && (is-il) < shortest) {
first = il;
last = is;
shortest = is - il;
}
if (is-il > longest) longest = is - il;
il = ++is;
}
if (longest > 1) {
m_shft (a, q, size, first, last);
m_hes (a, q, size);
}
}
return;
}