scilab-linear-algebra/jacobi.sce at master · saifsmailbox98/scilab-linear-algebra · GitHub

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

function [X, PN] = jacobi(A)
    if isequal(A,A.')

        [m, n] = size(A);
    check = %t;
    PN = eye(m, n);

    while check

        P = eye(m, n);
        p = 1;
        q = 2;
        for i=1:m
            for j=1:n
                if( i~=j & abs(A(i, j))>A(p,q))
                    p = i;
                    q = j;
                end
            end
        end


    t = 0.5 * atan(2*A(p, q)/(A(p, p) - A(q, q)));
    P(p, p) = cos(t);
    P(q, q) = cos(t);
    P(p, q) = -sin(t);
    P(q, p) = sin(t);

    A = P'*A*P;
    PN = PN * P;

    check = %f;
    for i=1:m
        for j=1:n
            if i~= j & abs(A(i, j) - 0) > 0.00000001
                check = %t;


end
end
end

end

X = diag(A)';
PN = PN';
end

endfunction


A = [5 0 1; 0 -2 0; 1 0 5];
[val, vec] = jacobi(A);

disp(val);
disp(vec);


A = [1 sqrt(2) 2; sqrt(2) 3 sqrt(2); 2 sqrt(2) 1];
[val, vec] = jacobi(A);

disp(val);
disp(vec);

/*

   6.  -2.   4.

   0.7071068   0.   0.7071068
   0.          1.   0.
  -0.7071068   0.   0.7071068

   5.   1.  -1.

   0.5         0.7071068   0.5
  -0.5         0.7071068  -0.5
  -0.7071068   0.          0.7071068

*/