---

from datetime import datetime
from pathlib import Path

import matplotlib.pyplot as plt
import numpy as np
from scipy.integrate import quad

import config

METHODS = (
    ("ЛПР", "Метод левых прямоугольников", "left_rectangles"),
    ("ППР", "Метод правых прямоугольников", "right_rectangles"),
    ("СПР", "Метод средних прямоугольников", "midpoint_rectangles"),
    ("ТР", "Метод трапеций", "trapezoids"),
    ("ММК", "Метод Монте-Карло для интеграла", "monte_carlo"),
    ("СИМП", "Метод Симпсона", "simpson"),
)


def f(x):
    """Вычисляем значение функции f(x) = sqrt((4-25*x**2)**3)."""
    x_values = np.asarray(x)
    return np.power(np.maximum(4 - 25 * x_values**2, 0), 1.5)


def true_area():
    """Возвращаем аналитическое (точное численное) значение площади."""
    val, _ = quad(f, config.A, config.B)
    return float(val)


def validate_config():
    """Проверяем корректность входных параметров."""
    if config.B <= config.A:
        raise ValueError("Параметр B должен быть больше A.")
    if not config.N_VALUES:
        raise ValueError("Параметр N_VALUES не должен быть пустым.")
    for n in config.N_VALUES:
        if int(n) <= 0:
            raise ValueError("Все значения N_VALUES должны быть больше нуля.")
    if int(config.PLOT_N) <= 0:
        raise ValueError("Параметр PLOT_N должен быть больше нуля.")


def make_rng(n):
    """Создаём генератор для метода Монте-Карло с устойчивым seed."""
    if getattr(config, "RANDOM_SEED", None) is None:
        return np.random.default_rng()
    return np.random.default_rng(int(config.RANDOM_SEED) + int(n))


def left_rectangles(n):
    """Вычисляем площадь методом левых прямоугольников."""
    dx = (config.B - config.A) / n
    x_left = config.A + np.arange(n) * dx
    return float(dx * np.sum(f(x_left)))


def right_rectangles(n):
    """Вычисляем площадь методом правых прямоугольников."""
    dx = (config.B - config.A) / n
    x_right = config.A + np.arange(1, n + 1) * dx
    return float(dx * np.sum(f(x_right)))


def midpoint_rectangles(n):
    """Вычисляем площадь методом средних прямоугольников."""
    dx = (config.B - config.A) / n
    x_mid = config.A + (np.arange(n) + 0.5) * dx
    return float(dx * np.sum(f(x_mid)))


def trapezoids(n):
    """Вычисляем площадь методом трапеций."""
    dx = (config.B - config.A) / n
    x_points = np.linspace(config.A, config.B, n + 1)
    y_points = f(x_points)
    return float(dx * (0.5 * y_points[0] + np.sum(y_points[1:-1]) + 0.5 * y_points[-1]))


def monte_carlo(n):
    """Вычисляем площадь методом Монте-Карло для интеграла."""
    rng = make_rng(n)
    x_random = rng.uniform(config.A, config.B, n)
    return float((config.B - config.A) * np.mean(f(x_random)))


def simpson(n):
    """Вычисляем площадь методом Симпсона для n параболических фигур."""
    subintervals = 2 * n
    dx = (config.B - config.A) / subintervals
    x_points = np.linspace(config.A, config.B, subintervals + 1)
    y_points = f(x_points)
    return float(
        dx
        / 3.0
        * (
            y_points[0]
            + y_points[-1]
            + 4.0 * np.sum(y_points[1:-1:2])
            + 2.0 * np.sum(y_points[2:-1:2])
        )
    )


METHOD_FUNCTIONS = {
    "left_rectangles": left_rectangles,
    "right_rectangles": right_rectangles,
    "midpoint_rectangles": midpoint_rectangles,
    "trapezoids": trapezoids,
    "monte_carlo": monte_carlo,
    "simpson": simpson,
}


def calculate_errors(estimate, exact_area):
    absolute_error = abs(estimate - exact_area)
    relative_error = absolute_error / exact_area * 100.0
    return absolute_error, relative_error


def build_results():
    exact_area = true_area()
    results = []
    for short_name, full_name, function_name in METHODS:
        row = {
            "short_name": short_name,
            "full_name": full_name,
            "values": {},
        }
        method = METHOD_FUNCTIONS[function_name]
        for n in config.N_VALUES:
            estimate = method(int(n))
            absolute_error, relative_error = calculate_errors(estimate, exact_area)
            row["values"][int(n)] = {
                "estimate": estimate,
                "absolute_error": absolute_error,
                "relative_error": relative_error,
            }
        results.append(row)
    return results


def format_number(value, digits=8, scientific_digits=10):
    """
    Форматирует число:
    - обычный fixed-point для стандартных значений,
    - scientific notation для очень маленьких чисел.
    """
    if value == 0:
        return f"{0:.{digits}f}"
    if abs(value) < 10 ** (-digits):
        return f"{value:.{scientific_digits}e}"
    return f"{value:.{digits}f}"


def format_result_cell(short_name, n, value):
    estimate_str = format_number(value["estimate"], digits=8, scientific_digits=10)

    if short_name in ("СПР", "ТР") and int(n) == 1000:
        absolute_str = format_number(value["absolute_error"], digits=16, scientific_digits=12)
        relative_str = format_number(value["relative_error"], digits=12, scientific_digits=12)
    else:
        absolute_str = format_number(value["absolute_error"], digits=8, scientific_digits=10)
        relative_str = format_number(value["relative_error"], digits=4, scientific_digits=10)

    return f"Ŝ={estimate_str}; Δ={absolute_str}; δ={relative_str}%"


def format_results_table(results):
    headers = ["Методы"] + [f"n={int(n)}" for n in config.N_VALUES]
    rows = []
    for row in results:
        rows.append(
            [row["short_name"]]
            + [format_result_cell(row["short_name"], int(n), row["values"][int(n)]) for n in config.N_VALUES]
        )

    widths = [
        max(len(headers[col]), *(len(row[col]) for row in rows))
        for col in range(len(headers))
    ]
    border = "+-" + "-+-".join("-" * width for width in widths) + "-+"

    lines = [border]
    lines.append(
        "| "
        + " | ".join(headers[col].ljust(widths[col]) for col in range(len(headers)))
        + " |"
    )
    lines.append(border)
    for row in rows:
        lines.append(
            "| "
            + " | ".join(row[col].ljust(widths[col]) for col in range(len(row)))
            + " |"
        )
    lines.append(border)
    return "\n".join(lines)


def find_best_results(results):
    best_by_n = {}
    overall_best = None

    for n in config.N_VALUES:
        n = int(n)
        best_row = min(
            results,
            key=lambda row: row["values"][n]["absolute_error"],
        )
        best_by_n[n] = {
            "short_name": best_row["short_name"],
            "full_name": best_row["full_name"],
            **best_row["values"][n],
        }

    for row in results:
        for n, value in row["values"].items():
            candidate = {
                "n": n,
                "short_name": row["short_name"],
                "full_name": row["full_name"],
                **value,
            }
            if (
                overall_best is None
                or candidate["absolute_error"] < overall_best["absolute_error"]
            ):
                overall_best = candidate

    return best_by_n, overall_best


def format_best_result(value):
    return (
        f"{value['short_name']} ({value['full_name']}): "
        f"Ŝ={format_number(value['estimate'], digits=10, scientific_digits=12)}; "
        f"Δ={format_number(value['absolute_error'], digits=16, scientific_digits=12)}; "
        f"δ={format_number(value['relative_error'], digits=12, scientific_digits=12)}%"
    )


def build_output_basename(output_dir):
    base = getattr(config, "PLOT_BASENAME", "plot")
    if getattr(config, "SAVE_UNIQUE_NAMES", False):
        timestamp = datetime.now().strftime(
            getattr(config, "FILENAME_TIMESTAMP_FORMAT", "%Y%m%d_%H%M%S")
        )
        base = f"{base}_{timestamp}"

    candidate = base
    index = 1
    formats = getattr(config, "SAVE_FORMATS", ["png"])
    while any((output_dir / f"{candidate}.{ext}").exists() for ext in formats):
        candidate = f"{base}_{index}"
        index += 1

    return candidate


def plot_function(ax):
    samples = getattr(config, "CURVE_SAMPLES", 400)
    x_line = np.linspace(config.A, config.B, samples)
    y_line = f(x_line)
    ax.plot(
        x_line,
        y_line,
        color=getattr(config, "FUNCTION_COLOR", "b"),
        linewidth=2.2,
        label=getattr(config, "FUNCTION_LABEL", r"$f(x) = \sqrt{(4-25x^2)^3}$"),
        zorder=3,
    )


def setup_method_axis(ax, title):
    ax.set_title(title, fontsize=10)
    ax.set_xlabel("x")
    ax.set_ylabel("y")
    ax.set_xlim(config.A, config.B)
    ax.set_ylim(0, float(np.max(f(np.linspace(config.A, config.B, 200)))) * 1.12)
    ax.grid(True, alpha=getattr(config, "GRID_ALPHA", 0.5))


def plot_left_rectangles(ax, n):
    dx = (config.B - config.A) / n
    x_left = config.A + np.arange(n) * dx
    ax.bar(
        x_left,
        f(x_left),
        width=dx,
        align="edge",
        color=getattr(config, "FIGURE_COLOR", "orange"),
        edgecolor=getattr(config, "FIGURE_EDGE_COLOR", "black"),
        alpha=0.35,
        linewidth=1.0,
        label="Левые прямоугольники",
    )


def plot_right_rectangles(ax, n):
    dx = (config.B - config.A) / n
    x_right = config.A + np.arange(1, n + 1) * dx
    ax.bar(
        x_right - dx,
        f(x_right),
        width=dx,
        align="edge",
        color=getattr(config, "FIGURE_COLOR", "green"),
        edgecolor=getattr(config, "FIGURE_EDGE_COLOR", "black"),
        alpha=0.35,
        linewidth=1.0,
        label="Правые прямоугольники",
    )


def plot_midpoint_rectangles(ax, n):
    dx = (config.B - config.A) / n
    x_left = config.A + np.arange(n) * dx
    x_mid = x_left + 0.5 * dx
    ax.bar(
        x_left,
        f(x_mid),
        width=dx,
        align="edge",
        color=getattr(config, "FIGURE_COLOR", "red"),
        edgecolor=getattr(config, "FIGURE_EDGE_COLOR", "black"),
        alpha=0.35,
        linewidth=1.0,
        label="Средние прямоугольники",
    )
    ax.scatter(
        x_mid,
        f(x_mid),
        s=18,
        color=getattr(config, "FIGURE_EDGE_COLOR", "black"),
        zorder=4,
    )


def plot_trapezoids(ax, n):
    x_points = np.linspace(config.A, config.B, n + 1)
    y_points = f(x_points)
    for i in range(n):
        polygon_x = [x_points[i], x_points[i], x_points[i + 1], x_points[i + 1]]
        polygon_y = [0, y_points[i], y_points[i + 1], 0]
        ax.fill(
            polygon_x,
            polygon_y,
            color=getattr(config, "FIGURE_COLOR", "purple"),
            edgecolor=getattr(config, "FIGURE_EDGE_COLOR", "black"),
            alpha=0.32,
            linewidth=1.0,
        )
    ax.plot(
        x_points,
        y_points,
        color=getattr(config, "FIGURE_EDGE_COLOR", "black"),
        linewidth=1.2,
        label="Трапеции",
    )


def plot_monte_carlo(ax, n):
    rng = make_rng(n)
    x_random = np.sort(rng.uniform(config.A, config.B, n))
    y_random = f(x_random)
    width = (config.B - config.A) / n
    ax.bar(
        x_random,
        y_random,
        width=width * 0.85,
        align="center",
        color=getattr(config, "MC_COLOR", "green"),
        edgecolor=getattr(config, "FIGURE_EDGE_COLOR", "black"),
        alpha=0.32,
        linewidth=0.9,
        label="Случайные прямоугольники",
    )
    ax.scatter(x_random, y_random, s=20, color=getattr(config, "MC_COLOR", "green"), zorder=4)


def plot_simpson(ax, n):
    subintervals = 2 * n
    x_points = np.linspace(config.A, config.B, subintervals + 1)
    y_points = f(x_points)
    for figure_index in range(n):
        i = 2 * figure_index
        local_x = x_points[i: i + 3]
        local_y = y_points[i: i + 3]
        coefficients = np.polyfit(local_x, local_y, deg=2)
        x_dense = np.linspace(local_x[0], local_x[-1], 80)
        y_dense = np.polyval(coefficients, x_dense)
        ax.fill_between(
            x_dense,
            0,
            y_dense,
            color=getattr(config, "SIMPSON_COLOR", "cyan"),
            alpha=0.28,
            edgecolor=getattr(config, "FIGURE_EDGE_COLOR", "black"),
            linewidth=0.8,
        )
        ax.plot(
            x_dense,
            y_dense,
            color=getattr(config, "SIMPSON_COLOR", "cyan"),
            linewidth=1.2,
        )
        ax.vlines(
            [local_x[0], local_x[-1]],
            0,
            [local_y[0], local_y[-1]],
            color=getattr(config, "FIGURE_EDGE_COLOR", "black"),
            alpha=0.55,
            linewidth=0.8,
        )
    ax.scatter(
        x_points,
        y_points,
        s=16,
        color=getattr(config, "FIGURE_EDGE_COLOR", "black"),
        zorder=4,
        label="Узлы",
    )


PLOTTERS = {
    "left_rectangles": plot_left_rectangles,
    "right_rectangles": plot_right_rectangles,
    "midpoint_rectangles": plot_midpoint_rectangles,
    "trapezoids": plot_trapezoids,
    "monte_carlo": plot_monte_carlo,
    "simpson": plot_simpson,
}


def plot_results():
    output_dir = Path(__file__).resolve().parent / getattr(config, "PLOT_DIR", "plots")
    output_dir.mkdir(parents=True, exist_ok=True)

    n = int(config.PLOT_N)
    fig, axes = plt.subplots(3, 2, figsize=getattr(config, "FIGURE_SIZE", (12, 10)))
    axes_flat = axes.ravel()
    exact_area = true_area()

    fig.suptitle(
        f"{getattr(config, 'PLOT_TITLE', 'Численное интегрирование')}\n"
        f"{getattr(config, 'FUNCTION_LABEL', 'f(x)')}, Sист={format_number(exact_area, digits=10, scientific_digits=12)}, n={n}",
        fontsize=13,
    )

    for ax, (short_name, full_name, function_name) in zip(axes_flat, METHODS):
        PLOTTERS[function_name](ax, n)
        plot_function(ax)
        estimate = METHOD_FUNCTIONS[function_name](n)
        absolute_error, relative_error = calculate_errors(estimate, exact_area)

        title_estimate = format_number(estimate, digits=8, scientific_digits=10)

        if short_name in ("СПР", "ТР") and n == 1000:
            title_abs = format_number(absolute_error, digits=16, scientific_digits=12)
            title_rel = format_number(relative_error, digits=12, scientific_digits=12)
        else:
            title_abs = format_number(absolute_error, digits=6, scientific_digits=10)
            title_rel = format_number(relative_error, digits=6, scientific_digits=10)

        setup_method_axis(
            ax,
            f"{short_name}: {full_name}\n"
            f"Ŝ={title_estimate}; Δ={title_abs}; δ={title_rel}%",
        )
        ax.legend(loc="best", fontsize=8)

    plt.tight_layout(rect=(0, 0, 1, 0.94))

    base_name = build_output_basename(output_dir)
    saved_paths = []
    formats = getattr(config, "SAVE_FORMATS", ["png"])
    for extension in formats:
        file_path = output_dir / f"{base_name}.{extension}"
        fig.savefig(file_path, dpi=getattr(config, "PLOT_DPI", 100), bbox_inches="tight")
        saved_paths.append(file_path)

    backend = plt.get_backend().lower()
    if getattr(config, "SHOW_PLOT", True) and "agg" not in backend:
        plt.show()
    else:
        plt.close(fig)

    return saved_paths


def print_report(results, saved_paths):
    exact_value = true_area()

    print("=" * 76)
    print("ЛАБОРАТОРНАЯ РАБОТА №4")
    print("Вычисление площади функции различными методами")
    print("=" * 76)
    print("f(x) = sqrt((4 - 25x^2)^3)")
    print(f"Отрезок интегрирования: [{config.A}, {config.B}]")
    print(f"Sист = {format_number(exact_value, digits=10, scientific_digits=12)}")
    print("-" * 76)
    print(format_results_table(results))
    print("-" * 76)
    best_by_n, overall_best = find_best_results(results)
    print("ЛУЧШИЙ МЕТОД ПО МИНИМАЛЬНОЙ АБСОЛЮТНОЙ ПОГРЕШНОСТИ")
    for n in config.N_VALUES:
        n = int(n)
        print(f"n={n}: {format_best_result(best_by_n[n])}")
    print("Лучший результат по всей таблице:")
    print(f"n={overall_best['n']}: {format_best_result(overall_best)}")
    print("=" * 76)
    print("Файлы графика:")
    for path in saved_paths:
        print(f"- {path}")


def main():
    try:
        validate_config()
    except ValueError as exc:
        raise SystemExit(f"Ошибка в параметрах config.py: {exc}") from exc

    results = build_results()
    saved_paths = plot_results()
    print_report(results, saved_paths)


if __name__ == "__main__":
    main()