"""Literal Python conversion of table3_3.m."""

import numpy as np


def table3_3(n_values=None):
    if n_values is None:
        n = np.array([4, 8, 12, 16, 160, 1600, 16000, 160000], dtype=int)
    else:
        n = np.asarray(n_values, dtype=int)
    N = len(n)
    print("\n n \t\t  Error", end="")
    for in_ in range(N):
        xe = np.ones(n[in_])
        a = 2 * np.ones(n[in_])
        b = np.ones(n[in_])
        c = np.ones(n[in_])
        f = 4 * np.ones(n[in_])
        f[0] = 3
        f[n[in_] - 1] = 3
        xc = tridiag(a, b, c, f)
        err = xe - xc
        eM = np.linalg.norm(err, np.inf)

        print("\n %d \t\t  %7.2e" % (n[in_], eM), end="")
    print("\n")


def tridiag(a, b, c, f):
    #  Solve the  N x N  tridiagonal system for y:
    #
    #  [ a(1)  c(1)                                  ] [  y(1)  ]   [  f(1)  ]
    #  [ b(2)  a(2)  c(2)                            ] [  y(2)  ]   [  f(2)  ]
    #  [       b(3)  a(3)  c(3)                      ] [        ]   [        ]
    #  [            ...   ...   ...                  ] [  ...   ] = [  ...   ]
    #  [                    ...    ...    ...        ] [        ]   [        ]
    #  [                        b(N-1) a(N-1) c(N-1) ] [ y(N-1) ]   [ f(N-1) ]
    #  [                                 b(N)  a(N)  ] [  y(N)  ]   [  f(N)  ]
    #
    #  f must be a vector (row or column) of length N
    #  a, b, c must be vectors of length N (note that b(1) and c(N) are not used)

    N = len(f)
    v = np.zeros(N)
    y = v.copy()

    w = a[0]
    y[0] = f[0] / w
    for i in range(1, N):
        v[i - 1] = c[i - 1] / w
        w = a[i] - b[i] * v[i - 1]
        y[i] = (f[i] - b[i] * y[i - 1]) / w

    for j in range(N - 2, -1, -1):
        y[j] = y[j] - v[j] * y[j + 1]

    return y


if __name__ == "__main__":
    table3_3()
