---

import numpy as np
import matplotlib.pyplot as plt

# Общие параметры
KAPPA = 0.5
X_MIN, X_MAX = -5.0, 5.0


def exact_solution(x, t, kappa=KAPPA):
    """Точное решение задачи Коши для u_t = kappa * u_xx."""
    return np.exp(-x**2 / (1 + 4 * kappa * t)) / np.sqrt(1 + 4 * kappa * t)


def solve_explicit(dx, T, kappa=KAPPA, r_target=0.4):
    """
    Явная разностная схема для уравнения теплопроводности.
    Возвращает: x, u, dt, число шагов по времени, фактическое r.
    """
    n_x = int(round((X_MAX - X_MIN) / dx)) + 1
    x = np.linspace(X_MIN, X_MAX, n_x)

    # Подбираем dt так, чтобы r было около r_target, а конечное время было ровно T
    dt0 = r_target * dx**2 / kappa
    n_t = int(np.ceil(T / dt0))
    dt = T / n_t
    r = kappa * dt / dx**2

    # Начальное условие u(x, 0) = exp(-x^2)
    u = exact_solution(x, 0.0, kappa)

    # Явная схема: обновляем только внутренние узлы
    for _ in range(n_t):
        u_new = np.zeros_like(u)
        u_new[1:-1] = u[1:-1] + r * (u[:-2] - 2 * u[1:-1] + u[2:])
        # Граничные условия: u = 0
        u = u_new

    return x, u, dt, n_t, r


def task1():
    print("\n=== Задание 1. Точное решение ===")
    x = np.linspace(X_MIN, X_MAX, 201)
    times = [0.0, 0.2, 0.5, 1.0]

    plt.figure(figsize=(9, 5))

    for t in times:
        u = exact_solution(x, t)
        umax = np.max(u)
        print(f"t = {t:.1f}: max u = {umax:.6f}")
        plt.plot(x, u, label=f"t = {t:.1f}")

    plt.title("Точное решение уравнения теплопроводности")
    plt.xlabel("Координата x")
    plt.ylabel("Температура u(x,t)")
    plt.grid(True)
    plt.legend()
    plt.tight_layout()
    plt.show()


def task2():
    print("\n=== Задание 2. Численное решение явной схемой ===")
    dx = 0.1
    T = 1.0

    x, u, dt, n_t, r = solve_explicit(dx, T, kappa=KAPPA, r_target=0.4)

    print(f"dx = {dx}")
    print(f"dt = {dt:.8f}")
    print(f"Число шагов по времени = {n_t}")
    print(f"Коэффициент r = {r:.6f}")
    print("Условие устойчивости r <= 0.5:", r <= 0.5)

    plt.figure(figsize=(9, 5))
    plt.plot(x, u, label="численное решение, T = 1.0")
    plt.title("Численное решение при T = 1.0")
    plt.xlabel("Координата x")
    plt.ylabel("Температура u")
    plt.grid(True)
    plt.legend()
    plt.tight_layout()
    plt.show()

    return x, u, T


def task3(x, u_num, T):
    print("\n=== Задание 3. Сравнение численного и точного решений ===")

    u_exact = exact_solution(x, T)
    error = np.abs(u_num - u_exact)
    max_error = np.max(error)

    print(f"Максимальная абсолютная ошибка = {max_error:.6e}")

    plt.figure(figsize=(9, 5))
    plt.plot(x, u_num, label="численное решение")
    plt.plot(x, u_exact, "--", label="точное решение")
    plt.title("Сравнение численного и точного решений при T = 1.0")
    plt.xlabel("Координата x")
    plt.ylabel("Температура u")
    plt.grid(True)
    plt.legend()
    plt.tight_layout()
    plt.show()

    plt.figure(figsize=(9, 5))
    plt.plot(x, error)
    plt.title("Абсолютная ошибка по координате")
    plt.xlabel("Координата x")
    plt.ylabel("|u_числ - u_точн|")
    plt.grid(True)
    plt.tight_layout()
    plt.show()


def task4():
    print("\n=== Задание 4. Влияние шага сетки ===")
    T = 1.0
    dx_values = [0.2, 0.1, 0.05]

    print(f"{'dx':>8} {'dt':>14} {'узлы':>8} {'шаги':>8} {'max error':>16} {'r':>10}")
    print("-" * 70)

    for dx in dx_values:
        x, u, dt, n_t, r = solve_explicit(dx, T, kappa=KAPPA, r_target=0.4)
        u_exact = exact_solution(x, T)
        max_error = np.max(np.abs(u - u_exact))

        print(
            f"{dx:8.2f} {dt:14.8f} {len(x):8d} {n_t:8d} "
            f"{max_error:16.6e} {r:10.6f}"
        )


def task5():
    print("\n=== Задание 5. Эксперимент с устойчивостью ===")
    dx = 0.1
    T = 0.5

    plt.figure(figsize=(9, 5))

    for name, r_target in [("A", 0.4), ("B", 0.6)]:
        x, u, dt, n_t, r = solve_explicit(dx, T, kappa=KAPPA, r_target=r_target)
        max_abs_u = np.max(np.abs(u))

        print(f"\nВариант {name}:")
        print(f"  r_target = {r_target}")
        print(f"  фактическое r = {r:.6f}")
        print(f"  dt = {dt:.8f}")
        print(f"  число шагов = {n_t}")
        print(f"  max|u| на последнем слое = {max_abs_u:.6e}")

        plt.plot(x, u, label=f"{name}: r = {r:.3f}, max|u| = {max_abs_u:.2e}")

    plt.title("Влияние устойчивости явной схемы при T = 0.5")
    plt.xlabel("Координата x")
    plt.ylabel("Температура u")
    plt.grid(True)
    plt.legend()
    plt.tight_layout()
    plt.show()


def main():
    task1()
    x, u_num, T = task2()
    task3(x, u_num, T)
    task4()
    task5()


if __name__ == "__main__":
    main()