import math


def f(x):
    #return math.exp(x)
    return math.sin(abs(x - 1/2))


def thomas_algorithm(n, k, y):
    a = [1] * n
    a[0] = 0
    b = [4] * n
    c = [1] * n
    c[n-1] = 0
    d = [0] * n
    for i in range(n):
        d[i] = k * (y[i+2] - 2 * y[i+1] + y[i])

    alpha = [0] * n
    beta = [0] * n
    x = [0] * n

    alpha[0] = - c[0] / b[0]
    beta[0] = d[0] / b[0]

    for i in range(1, n):
        alpha[i] = - c[i] / (a[i] * alpha[i-1] + b[i])
        beta[i] = (d[i] - a[i] * beta[i-1]) / (a[i] * alpha[i-1] + b[i])

    x[n-1] = beta[n-1]
    for i in range(n-2, -1, -1):
        x[i] = alpha[i] * x[i+1] + beta[i]

    return [0] + x + [0]


def coefficients_calculation(h, c, y, n):
    a = [0] * n
    b = [0] * n
    d = [0] * n
    for i in range(1, n):
        a[i-1] = y[i-1]
        b[i-1] = (y[i] - y[i-1]) / h - h * (c[i] + 2 * c[i-1]) / 3
        d[i-1] = (c[i] - c[i-1]) / (3 * h)

    return a, b, d


if __name__ == '__main__':
    x0 = 0
    xn = 1
    n = 10
    h = (xn - x0) / n

    x = [0] * (n+1)
    y = [0] * (n+1)

    for i in range(n+1):
        x[i] = x0 + i * h
        y[i] = f(x[i])

    c = thomas_algorithm(n-1, 3 / h ** 2, y)
    a, b, d = coefficients_calculation(h, c, y, n+1)

    with open('test.txt', 'w') as t:
        for i in range(2*n+1):
            xi = x0 + 0.5 * i * h
            ind = i // 2
            if ind == n:
                ind -= 1
            spline_x = a[ind] + b[ind] * (xi - x[ind]) + c[ind] * (xi - x[ind]) ** 2 + d[ind] * (xi - x[ind]) ** 3
            yi = f(xi)
            s = "x=%f | spl(x) = %f    |   f(x) = %f    |   %f" % (xi, spline_x, yi, spline_x - yi)
            t.write(s)
            t.write('\n')

    with open('test.txt','r') as t:
        for i, line in enumerate(t):
            print(i, line)