Advertisement
sashachca

Untitled

Apr 10th, 2018
173
0
Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
Octave 1.66 KB | None | 0 0
  1. function Eigen
  2.  
  3.   A = [2, -1, -8; -1, -3, 4; -8, 4, -3]
  4.   len = length(A);
  5.   [L, V] = jacobi(A, len);
  6.   eigenvalues = diag(L);
  7.   I = eye(len);
  8.  
  9.   printf('\nsobstvennye znachenia lambda:\n')
  10.   for i = 1:len
  11.     printf(' lambda%i = %i\n', i, eigenvalues(i));
  12.   end
  13.  
  14.   printf('\nsobstvennye vektory v:\n')
  15.   for i = 1:len
  16.     printf(' v%i = (%i, %i, %i)\n', i, V(1,i), V(2,i), V(3,i));
  17.   end
  18.  
  19.   printf('\nA*V = lambda*V \n');
  20.  
  21.   for i = 1:3
  22.     printf('\n %i-e:\n', i);
  23.     printf('\nA*v%i =\n', i);
  24.     disp(A*V(:,i));
  25.    
  26.     printf('\nlambda%i*v%i =\n', i, i);
  27.     disp(eigenvalues(i)*V(:,i));
  28.     printf('\n');
  29.   end
  30.  
  31.  
  32.  
  33.  
  34. end
  35.  
  36.  
  37. function [L, V] = jacobi(A, len)
  38.   t = 1;
  39.   e = 0.0001;
  40.   T = E = eye(len);
  41.  
  42.  
  43.   while t > e
  44.  
  45.     nondiagonal = A - diag(diag(A));
  46.     [a, I] = max(abs(nondiagonal));
  47.     [a, i] = max(a);
  48.     I = [I(i), i]; %нашли индексы максимального (по модулю) недиагонального элемента матрицы A
  49.     U = U_(A, I, len); % вычисляем матрицу U
  50.    
  51.     A = U'*A*U;
  52.     E *= U;
  53.    
  54.     t = t_(nondiagonal);
  55.   end
  56.   L = A;
  57.   V = E;
  58.  
  59. end
  60.  
  61.  
  62. function U = U_(A, I, len)
  63.   U = eye(len); % сначала U - единичная матрица размера A
  64.   i = I(1); j = I(2); % распаковываем индексы
  65.   fi = fi_(A, i, j); % считаем угол φ
  66.  
  67.   U(j, i) = sin(fi);
  68.   U(i, j) = -U(j, i);
  69.   U(i, i) = U(j, j) = cos(fi);
  70.  
  71. end
  72.  
  73.  
  74. function fi = fi_(A, i, j)
  75.   fi = 0;
  76.  
  77.   fi = 0.5*atan(2*A(i, j)/(A(i, i) - A(j, j)));
  78.  
  79. end
  80.  
  81.  
  82. function t = t_(ND)
  83.   t = 0;
  84.  
  85.   t = sqrt((sum(sum(ND.^2))));
  86.  
  87. end
Advertisement
Add Comment
Please, Sign In to add comment
Advertisement