---

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

# 1. Определение функции и пределов (Вариант 22, правый столбец)
def f(x):
    # np.maximum защищает от отрицательных значений под корнем из-за погрешностей float
    return np.power(np.maximum(4 - 25*x**2, 0), 1.5)

a, b = 0.2, 0.4
n_values = [10, 100, 1000]

# Истинное значение через scipy
true_val, _ = quad(f, a, b)

# --- Реализация методов ---

def left_rect(f, a, b, n):
    h = (b - a) / n
    return h * sum(f(a + i * h) for i in range(n))

def right_rect(f, a, b, n):
    h = (b - a) / n
    return h * sum(f(a + (i + 1) * h) for i in range(n))

def mid_rect(f, a, b, n):
    h = (b - a) / n
    return h * sum(f(a + h/2 + i * h) for i in range(n))

def trapezoid(f, a, b, n):
    h = (b - a) / n
    return h * (0.5 * (f(a) + f(b)) + sum(f(a + i * h) for i in range(1, n)))

def simpson(f, a, b, n):
    if n % 2 != 0: n += 1  # Для Симпсона n должно быть четным
    h = (b - a) / n
    x = np.linspace(a, b, n + 1)
    y = f(x)
    return h/3 * (y[0] + y[-1] + 4 * sum(y[1:-1:2]) + 2 * sum(y[2:-2:2]))

def monte_carlo(f, a, b, n):
    np.random.seed(42) # Для воспроизводимости
    x_rand = np.random.uniform(a, b, n)
    # Находим максимум функции на отрезке интегрирования
    f_max = np.max(f(np.linspace(a, b, 1000))) * 1.05 
    y_rand = np.random.uniform(0, f_max, n)
    under_curve = y_rand < f(x_rand)
    return (b - a) * f_max * (np.sum(under_curve) / n)

# --- Вычисления и вывод таблицы ---

methods = {
    "ЛПР": left_rect,
    "ППР": right_rect,
    "СПР": mid_rect,
    "ТР": trapezoid,
    "СИМП": simpson,
    "ММК": monte_carlo
}

print(f"{'n':<10} | {'Метод':<10} | {'Результат':<12} | {'Абс. Погрешность'}")
print("-" * 55)

for n in n_values:
    for name, func in methods.items():
        res = func(f, a, b, n)
        error = abs(true_val - res)
        print(f"{n:<10} | {name:<10} | {res:.8f} | {error:.8e}")
    print("-" * 55)

print(f"Истинное значение (True Area): {true_val:.8f}")

# --- Построение графика ---

x_plot = np.linspace(a - 0.05, b + 0.05, 400)
y_plot = f(x_plot)

# Данные для графика левых прямоугольников (как просили в задании n=10)
n_plot = 10
x_edges = np.linspace(a, b, n_plot + 1)
x_left = x_edges[:-1]
y_left = f(x_left)

plt.figure(figsize=(10, 6))
plt.plot(x_plot, y_plot, 'r', linewidth=2, label=r'$f(x) = \sqrt{(4 - 25x^2)^3}$')

# Фигуры ЛПР
plt.bar(x_left, y_left, width=(b-a)/n_plot, alpha=0.3, align='edge', 
        color='blue', edgecolor='black', label='Левые прямоугольники (n=10)')

plt.axhline(0, color='black', lw=1)
plt.axvline(0, color='black', lw=1)
plt.title('График функции и метод левых прямоугольников')
plt.xlabel('x')
plt.ylabel('y')
plt.legend()
plt.grid(True, linestyle='--', alpha=0.7)
plt.show()