---

import numpy as np
import matplotlib.pyplot as plt

# ============================================
# ВАРИАНТ 7
# ============================================

A = np.array([
    [7, -3, 12, -1],
    [3, 5, -6, 2],
    [2, 1, 5, 3],
    [-5, 9, 1, -11]
], dtype=float)

b = np.array([20, 11, 19, -33], dtype=float)

print("="*60)
print("ВАРИАНТ 7: Решение СЛАУ")
print("="*60)
print("\nМатрица A:")
print(A)
print("\nВектор b:", b)

# ============================================
# МЕТОД ЖОРДАНА-ГАУССА
# ============================================

def gauss_jordan(A, b):
    n = len(A)
    aug = np.hstack([A.copy(), b.reshape(-1,1)])
    
    print("\nНачальная расширенная матрица:")
    print(aug)
    
    for i in range(n):
        # Выбор главного элемента
        pivot_row = i
        max_val = abs(aug[i, i])
        for k in range(i+1, n):
            if abs(aug[k, i]) > max_val:
                max_val = abs(aug[k, i])
                pivot_row = k
        
        if pivot_row != i:
            aug[[i, pivot_row]] = aug[[pivot_row, i]]
            print(f"\nШаг {i+1}: перестановка строк {i+1} и {pivot_row+1}")
        
        aug[i] = aug[i] / aug[i, i]
        print(f"\nШаг {i+1}.1: нормализация строки {i+1}")
        
        for j in range(n):
            if j != i:
                aug[j] = aug[j] - aug[j, i] * aug[i]
        print(f"Шаг {i+1}.2: исключение из остальных строк")
        print(aug)
    
    return aug[:, -1]

solution_gj = gauss_jordan(A, b)

print("\n" + "="*40)
print("РЕШЕНИЕ МЕТОДОМ ЖОРДАНА-ГАУССА:")
print("="*40)
for i, val in enumerate(solution_gj, 1):
    print(f"x{i} = {val:.10f}")

# Проверка
residual = A @ solution_gj - b
print("\nНевязки:")
for i, r in enumerate(residual, 1):
    print(f"Уравнение {i}: {r:.2e}")

# ============================================
# МЕТОД ЗЕЙДЕЛЯ
# ============================================

print("\n" + "="*60)
print("МЕТОД ЗЕЙДЕЛЯ")
print("="*60)

def seidel_method(A, b, x0=None, tol=1e-12, max_iter=100):
    n = len(A)
    if x0 is None:
        x0 = np.zeros(n)
    
    x = x0.copy()
    history = [x.copy()]
    
    print("\n Итерация         x1              x2              x3              x4          погрешность")
    print("-" * 85)
    
    for k in range(max_iter):
        x_prev = x.copy()
        
        for i in range(n):
            s = 0
            for j in range(n):
                if j != i:
                    s += A[i, j] * x[j]
            x[i] = (b[i] - s) / A[i, i]
        
        err = np.linalg.norm(x - x_prev, np.inf)
        history.append(x.copy())
        
        print(f"{k+1:3d}     {x[0]:12.10f}   {x[1]:12.10f}   {x[2]:12.10f}   {x[3]:12.10f}   {err:.2e}")
        
        if err < tol:
            print(f"\nСходимость достигнута за {k+1} итераций")
            break
    
    return x, history, k+1

zeidel_solution, history, iterations = seidel_method(A, b, tol=1e-12, max_iter=50)

print("\n" + "="*40)
print("СРАВНЕНИЕ РЕЗУЛЬТАТОВ:")
print("="*40)
print("Жордан-Гаусс:", [f"{x:.10f}" for x in solution_gj])
print("Зейдель:    ", [f"{x:.10f}" for x in zeidel_solution])

# ============================================
# МАТРИЦА ИТЕРАЦИЙ ЗЕЙДЕЛЯ
# ============================================

print("\n" + "="*60)
print("МАТРИЦА ИТЕРАЦИЙ МЕТОДА ЗЕЙДЕЛЯ")
print("="*60)

def seidel_iteration_matrix(A):
    n = len(A)
    L = np.zeros((n, n))
    D = np.zeros((n, n))
    R = np.zeros((n, n))
    
    for i in range(n):
        for j in range(n):
            if i > j:
                L[i, j] = A[i, j]
            elif i == j:
                D[i, j] = A[i, j]
            else:
                R[i, j] = A[i, j]
    
    print("\nL (нижняя треугольная):")
    print(L)
    print("\nD (диагональная):")
    print(D)
    print("\nR (верхняя треугольная):")
    print(R)
    
    LD = L + D
    print("\n(L + D):")
    print(LD)
    
    LD_inv = np.linalg.inv(LD)
    print("\n(L + D)^(-1):")
    print(LD_inv)
    
    B_seidel = -np.dot(LD_inv, R)
    print("\n" + "="*40)
    print("ИТЕРАЦИОННАЯ МАТРИЦА ЗЕЙДЕЛЯ B_seidel:")
    print("="*40)
    print(B_seidel)
    
    norm_B = np.max(np.sum(np.abs(B_seidel), axis=1))
    print(f"\nНорма матрицы Зейделя: ||B|| = {norm_B:.6f}")
    
    if norm_B < 1:
        print(f"\n✓ Норма < 1 => метод Зейделя СХОДИТСЯ")
        print(f"  Оценка скорости: 1/||B|| = {1/norm_B:.2f}")
    else:
        print(f"\n✗ Норма >= 1 => метод может расходиться")
    
    return B_seidel, norm_B

B_seidel, norm_B = seidel_iteration_matrix(A)

# ============================================
# ГРАФИК СХОДИМОСТИ (ДОЛЖЕН РАБОТАТЬ)
# ============================================

print("\n" + "="*60)
print("ПОСТРОЕНИЕ ГРАФИКА СХОДИМОСТИ...")
print("="*60)

# Вычисляем погрешности на каждой итерации
errors = []
for k in range(len(history)):
    error = np.linalg.norm(history[k] - solution_gj, np.inf)
    errors.append(error)

# Создаём график
plt.figure(figsize=(10, 6))
plt.semilogy(errors, 'b-o', linewidth=2, markersize=6)
plt.xlabel('Номер итерации', fontsize=12)
plt.ylabel('Погрешность (логарифмическая шкала)', fontsize=12)
plt.title('Сходимость метода Зейделя (Вариант 7)', fontsize=14)
plt.grid(True, alpha=0.3)
plt.axhline(y=1e-12, color='r', linestyle='--', label='Точность 1e-12')
plt.legend()
plt.tight_layout()

# Сохраняем и показываем
plt.savefig('convergence_variant7.png', dpi=150)
plt.show()

print("\nГрафик сохранен: convergence_variant7.png")
print("График должен открыться в отдельном окне!")

# ============================================
# ДОПОЛНИТЕЛЬНЫЙ ВАРИАНТ (по желанию)
# ============================================

print("\n" + "="*60)
print("РАБОТА ЗАВЕРШЕНА")
print("="*60)