"""Literal Python conversion of fig10_19.m."""

from __future__ import annotations

import matplotlib.pyplot as plt
import numpy as np

from subaxis import subaxis


def _set_bold_tick_font(ax: plt.Axes, fontsize: int) -> None:
    ax.tick_params(labelsize=fontsize)
    ax.xaxis.label.set_size(fontsize)
    ax.yaxis.label.set_size(fontsize)
    ax.xaxis.label.set_fontweight("bold")
    ax.yaxis.label.set_fontweight("bold")
    for label in ax.get_xticklabels() + ax.get_yticklabels():
        label.set_fontweight("bold")


def _print_latex_matrix(A: np.ndarray) -> None:
    for row_index, row in enumerate(A):
        ending = r" \\[0.2em]" if row_index < A.shape[0] - 1 else r" \\"
        print(
            " "
            + r" \,\,&\,\, ".join(f"{value:8.2f}" for value in row)
            + ending
        )
    print()


def fig10_19() -> None:
    # random generation of A values
    # m=10;   # number of vertices (must be even)
    # A[0:m//2,0:m//2]=10+np.random.randint(1,21,(m//2,m//2))
    # A[m//2:m,m//2:m]=10+np.random.randint(1,21,(m//2,m//2))
    # A[m//2,0]=1
    # A=A/np.sum(A,axis=0)

    # book example
    m = 4
    A = np.zeros((m, m))
    d = 1 / 20
    A[0, 0] = 1 / 4 - d
    A[1, 0] = 3 / 4 - d
    A[2, 0] = 2 * d
    A[0, 1] = 1 / 3
    A[1, 1] = 2 / 3
    A[2, 2] = 1 / 4
    A[3, 2] = 3 / 4
    A[2, 3] = 1 / 3
    A[3, 3] = 2 / 3

    # np.sum(A,axis=0)
    print("A =\n", A)
    eigens = np.linalg.eigvals(A)
    print("eigens =", eigens)
    sings = np.linalg.svd(A, compute_uv=False)
    print("sings =", sings)

    T = 40
    k = 0
    N = np.zeros(T - m, dtype=int)
    error = np.zeros(T - m)
    for im in range(m + 1, T + 1):
        # for im in range(4,T+1):
        k = k + 1
        N[k - 1] = im  # number of time steps
        y = np.zeros((N[k - 1], m))
        y[0, 0] = 100
        # y[0,:]=np.ones(m)
        for i in range(1, N[k - 1]):
            z = A @ y[i - 1, :].T
            y[i, :] = z.T
        q = 0.1 * (2 * np.random.rand(m) - 1)
        # q=0
        yd = (1 + q) * y
        XD = yd[0 : N[k - 1] - 1, :]
        YD = yd[1 : N[k - 1], :]
        n = N[k - 1] - 1

        # use the pseudo-inverse
        # X=XD
        # Y=YD
        # G=X.T@X
        # GI=np.linalg.inv(G)
        # kappa[k-1]=np.linalg.cond(G,np.inf)
        # P=Y.T@X@GI

        # use the SVD
        U, singular_values, Vh = np.linalg.svd(XD, full_matrices=True)
        S = np.zeros((n, m))
        S[:m, :m] = np.diag(singular_values)
        V = Vh.T
        SS = np.zeros((n, m))
        for i in range(m):
            SS[i, i] = 1 / S[i, i]
        P = YD.T @ U @ SS @ V.T

        Ac = P
        error[k - 1] = np.linalg.norm(A - Ac, np.inf) / np.linalg.norm(A, np.inf)

        if N[k - 1] == 5:
            # eigenvalues, VV=np.linalg.eig(P)
            # DD=np.diag(eigenvalues)
            # VI=np.linalg.inv(VV)
            # DD[3,3]*np.outer(VV[:,3],VI[3,:])+DD[2,2]*np.outer(VV[:,2],VI[2,:])+DD[1,1]*np.outer(VV[:,1],VI[1,:])+DD[0,0]*np.outer(VV[:,0],VI[0,:])-P
            # np.outer(VV[:,3],VI[3,:])
            # np.outer(VV[:,2],VI[2,:])
            # np.outer(VV[:,1],VI[1,:])
            # np.outer(VV[:,0],VI[0,:])

            _print_latex_matrix(A)
            _print_latex_matrix(Ac)

    print("max(error) =", np.max(error))

    fig1 = plt.figure(1, figsize=(7.2, 3.4))
    fig1.clf()
    # MATLAB figure position was [25 1115 658 230]; figsize is scaled to fit the plot window.
    ax1 = subaxis(1, 1, 1, 1, "MT", 0.01, "MB", 0.18, "MR", -0.02, "ML", 0.05, "P", 0.04)
    co = np.array(
        [
            [0, 0, 1],
            [0, 0.5, 0],
            [1, 0, 0],
            [0, 0.75, 0.75],
            [0.75, 0, 0.75],
            [0.75, 0.75, 0],
            [0.25, 0.25, 0.25],
        ]
    )
    ax1.set_prop_cycle(color=co.tolist())
    ax1.plot(N, error, "--r", linewidth=0.8)
    ax1.plot(N, error, "or", linewidth=1.5, markersize=8)
    ax1.set_frame_on(True)
    ax1.grid(True)
    # ax1.axis([a,b,-2.5,4])
    ax1.set_yticks([0, 0.05, 0.1])
    ax1.set_xlabel("n-axis", fontweight="bold")
    ax1.set_ylabel("Relative Error", fontweight="bold")
    _set_bold_tick_font(ax1, 14)  # MATLAB FontSize 16 reduced by 2

    # fig1.savefig('carserror.eps',format='eps',bbox_inches='tight')

    # plot y values at location m
    fig2 = plt.figure(2, figsize=(7.2, 3.4))
    fig2.clf()
    # MATLAB figure position was [27 797 658 230]; figsize is scaled to fit the plot window.
    ax2 = subaxis(1, 1, 1, 1, "MT", 0.01, "MB", 0.18, "MR", -0.02, "ML", 0.06, "P", 0.04)
    ax2.set_prop_cycle(color=co.tolist())
    ax2.plot(np.arange(1, N[k - 1]), YD[:, m - 1], "--ob")

    plt.show()


if __name__ == "__main__":
    fig10_19()
