№9 Разностный метод решения краевой задачи для ОДУ 2-го порядка

import math
def f(x):
    V = 8
    return 4 * pow(x, 4) - 3 * V * pow(x, 3) + 6 * x - 2 * V

def p(x):
    return x * x

def q(x):
    return x

def check_y(x):
    V = 8
    return x * x * (x - V)

def tridiagonal_algorithm(a, b, c, d):
    n = len(a)
    p = [0] * (n + 1)
    q = [0] * (n + 1)
    for i in range(n):
        p[i + 1] = c[i] / (b[i] - a[i] * p[i])
        q[i + 1] = (a[i] * q[i] - d[i]) / (b[i] - a[i] * p[i])

    res = [0] * n
    res[n - 1] = q[n]
    for i in range(n - 2, -1, -1):
        res[i] = p[i + 1] * res[i + 1] + q[i + 1]

    return res

def print_result(vx, vy):
    k = 10
    for i in range(0, len(vx), k):
        m = min(i + k, len(vx))
        for j in range(i, m):
            print(f"{vx[j]:.5f}", end="\t")
        print()
        for j in range(i, m):
            print(f"{vy[j]:.5f}", end="\t")
        print()
        for j in range(i, m):
            print(f"{check_y(vx[j]):.5f}", end="\t")
        print()
        for j in range(i, m):
            print(f"{vy[j] - check_y(vx[j]):.5f}", end="\t")
        print()
        print()

    mx = 0
    for i in range(len(vy)):
        mx = max(mx, abs(vy[i] - check_y(vx[i])))
    print(f"Максимальное отклонение: {mx:.5f}")

def solve():
    n = 6
    l, r = 0, 16  # T = V = 16
    h = (r - l) / n

    vx = [l + i * h for i in range(n + 1)]

    a = [0] * (n - 1)
    b = [0] * (n - 1)
    c = [0] * (n - 1)
    d = [0] * (n - 1)

    for i in range(n - 1):
        x = vx[i + 1]
        a[i] = (1 / h / h - p(x) / 2 / h)
        b[i] = -(-2 / h / h + q(x))
        c[i] = (1 / h / h + p(x) / 2 / h)
        d[i] = f(x)

    a[0] = 0
    c[-1] = 0

    vy = tridiagonal_algorithm(a, b, c, d)
    vy.insert(0, 0)
    vy.append(0)

    print_result(vx, vy)


if __name__ == "__main__":
    solve()