import numpy as np
import matplotlib.pyplot as plt
from cycler import cycler

from plot_helpers import CO, bold_ticks


def fig7_13():
    # using velocity Verlet on the system y' = f(t,y) with y(0) = y0

    t0 = 0
    y0 = np.array([np.pi / 4, 0], dtype=float)
    tmax = 4000
    k = 1 / 2
    n = round(tmax / k)
    print(n)
    t = np.linspace(t0, tmax, n)

    y, v = verlet(t, y0, k, n)

    # determine the number of zero crossings
    # k=1/2 => 1009 crossing for 0<t<3266
    # nzeros=0
    # for it in range(1,n):
    #     if y[it]*y[it-1] < 0:
    #         nzeros=nzeros+1
    #
    # nzeros_verlet=nzeros

    fig1, ax1 = plt.subplots(figsize=(8.8, 3.3))
    fig1.subplots_adjust(left=0.09, right=0.98, bottom=0.22, top=0.98)
    ax1.set_prop_cycle(cycler(color=CO))

    ax1.plot(t, y, "-b", linewidth=2)
    T1 = 40
    ax1.set_xlim(0, T1)
    ax1.set_ylim(-0.85, 0.85)
    ax1.plot([0, T1], [np.pi / 4, np.pi / 4], "--r", linewidth=2)
    ax1.set_yticks([-np.pi / 4, 0, np.pi / 4])
    ax1.set_yticklabels([" ", "0", " "])
    say = r"$-\pi/4$"
    ax1.text(-3, -np.pi / 4, say, fontsize=18)
    say = r"$\pi/4$"
    ax1.text(-2.6, np.pi / 4, say, fontsize=18)
    ax1.set_xticks([0, 10, 20, 30, 40])
    ax1.grid(True)
    ax1.set_xlabel("t-axis")
    ax1.set_ylabel("Theta")
    bold_ticks(ax1, 14)
    for spine in ax1.spines.values():
        spine.set_visible(True)

    # plt.savefig('/Users/mark/Desktop/ver1.eps')

    fig2, ax2 = plt.subplots(figsize=(8.8, 3.3))
    fig2.subplots_adjust(left=0.08, right=0.995, bottom=0.22, top=0.98)
    ax2.set_prop_cycle(cycler(color=CO))

    tr = np.linspace(0, tmax, 100)
    aa = 1e-4 * np.log(7 * np.pi / 3)
    R = np.zeros(100)
    for ir in range(100):
        R[ir] = 0.25 * np.pi * np.exp(-aa * tr[ir])

    ax2.plot(t, y, "-b", linewidth=0.01, color="#00BFFF")
    ax2.plot([0, tmax], [np.pi / 4, np.pi / 4], "--r", linewidth=2)
    ax2.set_yticks([-np.pi / 4, 0, np.pi / 4])
    # ax2.set_yticklabels(["-pi/4", "0", "pi/4"])
    ax2.set_yticklabels([" ", "0", " "])
    say = r"$-\pi/4$"
    ax2.text(-300, -np.pi / 4, say, fontsize=18)
    say = r"$\pi/4$"
    ax2.text(-250, np.pi / 4, say, fontsize=18)
    ax2.grid(True)
    ax2.set_xlabel("t-axis")
    ax2.set_ylabel("Theta")
    bold_ticks(ax2, 14)
    for spine in ax2.spines.values():
        spine.set_visible(True)

    # plt.savefig('/Users/mark/Desktop/ver2.eps')
    plt.show()


# right-hand side of DE
def f(t, y):
    z = np.array([y[1], -np.sin(y[0])])
    return z


# Verlet method
def verlet(t, y0, k, n):
    ypoints = np.empty(n)
    vpoints = np.empty(n)
    ypoints[0] = y0[0]
    vpoints[0] = y0[1]
    y = y0[0]
    v = y0[1]
    for i in range(1, n):
        Y = y + k * v - 0.5 * k ** 2 * np.sin(y)
        V = v - 0.5 * k * np.sin(y) - 0.5 * k * np.sin(Y)
        ypoints[i] = Y
        vpoints[i] = V
        y = Y
        v = V
    return ypoints, vpoints


if __name__ == "__main__":
    fig7_13()
