Метод вращения. Метод Якоби. Поиск собственных значений матрицы.

[@OnTheWayToTheDirection]

Команды для генерации анимации:
manim main.py play_whole_scenario - полный рендер
manim main.py -pql - быстрый рендер

[@arsenozd]

А что делать, если диагональные элементы матрицы равны дур другу?

[@OnTheWayToTheDirection]

Python Manim скрипт, с помощью которого созданы некоторые анимации:
from manim import *
import math
from math import cos
from math import sin
class MethodYacoby(MovingCameraScene):
def construct(self):
A = np.array([[1, 2, -1],
[2, 3, 4],
[-1, 4, 5]])

eps = 0.0001 # точность
n = 3
eigenvalues = np.zeros(n) # Инициализируем вектор собственных значений


#-
self.play(
Write(Text("Матрицы А:").to_corner(UL).scale(0.5)),
Write(Text("Матрицы H:").to_corner(UR).scale(0.5))
)
#-
# Повторяем вращения Якоби, пока не достигнем требуемой точности
iteration = 0
while True:
# Находим наибольший по модулю внедиагональный элемент
max_off_diag = 0
p = q = 0
for i in range(n):
for j in range(i+1, n):
if abs(A[i, j]) > max_off_diag:
max_off_diag = abs(A[i, j])
p, q = i, j

#-
_A = Matrix(A.round(2)).scale(0.6).move_to(LEFT*4.7 + DOWN*3*iteration)
matrix_label = Text(f"Матрица A:",font=BOLD, font_size=13).next_to(_A, UP)
i = Text(f"Итерация: {iteration}",font=BOLD, font_size=13).next_to(_A, LEFT).rotate(PI/2)
self.play(Write(i))
self.play(Write(matrix_label))
self.play(Write(_A))
#-
# Если внедиагональные элементы пренебрежимо малы, мы достигли сходимости
if max_off_diag < eps:
eigenvalues = np.diag(A)
break
# Вычисляем угол вращения Якоби
theta = 0.5 * math.atan2(2 * A[p, q], A[p, p] - A[q, q])
#-
inf_label = Text(f"Индексы максимального недиагонального элемента

p = max_i = {p}

q = max_j = {q}

Максимальный недиагональный элемент: {max_off_diag}", font=BOLD, font_size=13).move_to(DOWN*3*iteration + UP/2)
self.play(Write(inf_label))
#-
#-
ang_label = Text(f"Угол поворота = {round(theta, 2)} [радиан]",font=BOLD, font_size=13).next_to(inf_label, DOWN)
self.play(Write(ang_label))
#-
#-
error_label = Text(f"Погрешность: {max_off_diag}, Требуемая точность: {eps}",font=BOLD, font_size=13).next_to(ang_label, DOWN)
self.play(Write(error_label))
#-
# Строим матрицу вращения H
H = np.eye(n) # ?
"""
H = np.eye(3) = [[1, 0, 0],
[0, 1, 0],
[0, 0, 1]]
"""
H[p, p] = H[q, q] = math.cos(theta)
H[p, q] = -math.sin(theta)
H[q, p] = math.sin(theta)
#-
_H = Matrix(H.round(2)).scale(0.65).move_to(RIGHT*5.5 + DOWN*3*iteration)
matrix_label = Text(f"Матрица поворота H:",font=BOLD, font_size=13).next_to(_H, UP)
self.play(Write(matrix_label))
self.play(Write(_H))
#-
self.play(
self.camera.frame.animate.move_to(DOWN*3*(iteration+1))
)
# Применяем вращение к матрице A
A = H.T @ A @ H
iteration += 1
eigenvalues[0]
eigenvalues[1]
eigenvalues[2]

[@OnTheWayToTheDirection]

Алгоритм Якоби python:
```python
import numpy as np
import math
def jacobi_eigenvalues_H(A, eps=1e-10):
n = 3
eigenvalues = np.zeros(n) # Инициализируем вектор собственных значений
iteration = 0
# Повторяем вращения Якоби, пока не достигнем требуемой точности
while True:
# Находим наибольший по модулю внедиагональный элемент
max_off_diag = 0
p = q = 0
for i in range(n):
for j in range(i+1, n):
if abs(A[i, j]) > max_off_diag:
max_off_diag = abs(A[i, j])
p, q = i, j
# Если внедиагональные элементы пренебрежимо малы, мы достигли сходимости
if max_off_diag < eps:
eigenvalues = np.diag(A)
break
# Вычисляем угол вращения Якоби
theta = 0.5 * math.atan2(2 * A[p, q], A[p, p] - A[q, q])
# Строим матрицу вращения H
H = np.eye(n) # ?
H[p, p] = H[q, q] = math.cos(theta)
H[p, q] = -math.sin(theta)
H[q, p] = math.sin(theta)
print("Матрица поворота:")
print(H)
# Применяем вращение к матрице A
A = H.T @ A @ H
iteration += 1
return eigenvalues
A = np.array([[1, 2, -1],
[2, 3, 4],
[-1, 4, 5]])
eigenvalues = jacobi_eigenvalues_H(A)
print("Собственные значения:", eigenvalues)
```