Загрузка данных


import sympy as sp

x = sp.symbols("x")
f_expr = 2 * sp.log(x) - x ** 2 + 5                      # уравнение f(x) = 0

SEGMENTS = [(0.07, 0.09), (2, 3)]                        # отрезки изоляции корней (по графику)
ROOTS_COUNT = 2                                          # сколько корней найти
EPSILON = 0.0001
PARTS = 1000                                             # на сколько частей делим отрезок (для min и max)

if ROOTS_COUNT > len(SEGMENTS):
    raise ValueError(f"В уравнении корней: {len(SEGMENTS)}, а запрошено {ROOTS_COUNT}")

f = sp.lambdify(x, f_expr, "math")
f1 = sp.lambdify(x, sp.diff(f_expr, x), "math")
f2 = sp.lambdify(x, sp.diff(f_expr, x, 2), "math")


def get_points(A, B):
    return [A + (B - A) * k / PARTS for k in range(PARTS + 1)]


def min_max(g, A, B):                                    # min и max |g(x)| на [A, B]
    values = [abs(g(p)) for p in get_points(A, B)]
    return min(values), max(values)


def constant_sign(g, A, B):                              # g не меняет знак на [A, B]
    values = [g(p) for p in get_points(A, B)]
    return all(v > 0 for v in values) or all(v < 0 for v in values)


def start_point(A, B):                                   # начальный конец: где f(x) * f''(x) > 0
    if f(A) * f2(A) > 0:
        return A, B
    return B, A


def dichotomy(A, B):
    left, right = A, B
    value_left = f(left)

    if value_left * f(right) > 0:                        # проверка знаков до цикла
        return None

    previous_x = None
    iterations = 0
    while True:
        x_new = (left + right) / 2
        value_x = f(x_new)
        iterations += 1

        if value_x == 0:                                 # точное попадание в корень
            break
        if previous_x is not None and abs(x_new - previous_x) < EPSILON:   # |xn - xn-1| < eps
            break

        if value_left * value_x < 0:
            right = x_new
        else:
            left = x_new
            value_left = value_x
        previous_x = x_new

    return x_new, iterations


def simple_iterations(A, B):
    lam = -1 / f1((A + B) / 2)                           # lambda = -1 / f'(x0), x0 - середина отрезка
    phi = lambda t: t + lam * f(t)                       # x = phi(x) = x + lambda * f(x)
    q = min_max(lambda t: 1 + lam * f1(t), A, B)[1]      # q = max |phi'(x)|

    if q >= 1:                                           # условие сходимости не выполнено
        return None

    x_prev = A
    iterations = 0
    while True:
        x_new = phi(x_prev)
        iterations += 1
        if q / (1 - q) * abs(x_new - x_prev) < EPSILON:  # q/(1-q) * |xn - xn-1| < eps
            break
        x_prev = x_new

    return x_new, iterations


def newton(A, B):
    if not (constant_sign(f1, A, B) and constant_sign(f2, A, B) and f1(A) * f2(A) > 0):
        return None                                      # f' и f'' должны иметь одинаковый знак

    m1 = min_max(f1, A, B)[0]                            # m = min |f'(x)|
    x_prev = start_point(A, B)[0]
    iterations = 0
    while True:
        x_new = x_prev - f(x_prev) / f1(x_prev)
        iterations += 1
        if abs(f(x_new)) / m1 < EPSILON:                 # |f(xn)| / m < eps
            break
        x_prev = x_new

    return x_new, iterations


def secant(A, B):
    if not (constant_sign(f1, A, B) and constant_sign(f2, A, B)):
        return None

    m1 = min_max(f1, A, B)[0]
    start, other = start_point(A, B)
    x_prev, x_cur = other, start                         # нужны две начальные точки
    iterations = 0
    while True:
        x_new = x_cur - f(x_cur) * (x_cur - x_prev) / (f(x_cur) - f(x_prev))
        iterations += 1
        if abs(f(x_new)) / m1 < EPSILON:                 # критерий как у Ньютона
            break
        x_prev, x_cur = x_cur, x_new

    return x_new, iterations


def chord(A, B):
    if not (constant_sign(f1, A, B) and constant_sign(f2, A, B)):
        return None

    m1 = min_max(f1, A, B)[0]
    fixed, x_cur = start_point(A, B)                     # fixed - неподвижный конец хорды
    iterations = 0
    while True:
        x_new = x_cur - f(x_cur) * (fixed - x_cur) / (f(fixed) - f(x_cur))
        iterations += 1
        if abs(f(x_new)) / m1 < EPSILON:                 # критерий как у Ньютона
            break
        x_cur = x_new

    return x_new, iterations


def show(name, result):
    if result is None:
        print(f"  {name}: метод неприменим на этом отрезке")
    else:
        print(f"  {name}: {result[0]:.8f}, итераций: {result[1]}")


for i in range(ROOTS_COUNT):
    A, B = SEGMENTS[i]
    print(f"Корень {i + 1}, отрезок [{A}; {B}]")
    show("Дихотомия       ", dichotomy(A, B))
    show("Простые итерации", simple_iterations(A, B))
    show("Ньютон          ", newton(A, B))
    show("Секущие         ", secant(A, B))
    show("Хорды           ", chord(A, B))