import matplotlib.pyplot as plt
import numpy as np
import scipy.integrate as si
import scipy.interpolate as interp
from array import *
import math

#f(x,y)
def f(x,y):
    return -y - x**2

def Adams(x, y, num, h):
    for i in range(num-1):
        y.insert(len(y), Euler(y,x,h))
    F = [f(x[0],y[0]), f(x[1],y[1]), f(x[2],y[2]), f(x[3],y[3])]
    coef = 24
    for i  in range(len(x)):
        if i >=  (len(y)):
            y.insert(i, y[i-1] + h/coef*(55*F[3] - 59*F[2] + 37*F[1] - 9*F[0]))
            for j in range(3):
                F[j] = F[j+1]
            F[3] = f(x[i],y[i])
    return y


#Метод эйлера
def Euler(y, x, h):
  return y[len(y)-1] + h*f(x[len(y) - 1], y[len(y) - 1]) 

     
def RRKM(b, h, e):
     y = array('f', [10])
     x = array('f', [0])
     k = array('f', [0,0,0,0,0,0])
     flag = 0

     while x[len(x) - 1]<=b:
        while True:
            k[1] = h*f(x[len(y) - 1], y[len(y) - 1])
            k[2] = h*f(x[len(y) - 1] + 1/3*h, y[len(y) - 1] +1/3*k[1])
            k[3] = h*f(x[len(y) - 1] + 1/3*h, y[len(y) - 1] + 1/6*k[1] + 1/6*k[2])
            k[4] = h*f(x[len(y) - 1] + 1/2*h, y[len(y) - 1] + 1/8*k[1] + 3/8*k[3])
            k[5] = h*f(x[len(y) - 1] + h, y[len(y) - 1] + 1/2*k[1] - 3/2*k[3] + 2*k[4])
            q=(8*k[4] - 9*k[3] - k[5] + 2*k[1])/30
            if abs(q)<e*abs(y[len(y) - 1]):
                break
            h = h/2
        x.insert(len(x), x[len(x) - 1] + h)
        y.insert(len(y), y[len(y) - 1] + 1/6*k[1] + 2/3*k[4] + 1/6*k[5])
        if abs(q)<=e*abs(y[len(y) - 1])/32:
            h = h*2
     return x,y

def RK4(x, h, b, n):
    y = array('f', [10])
    k = array('f', [0,0,0,0,0])
    while len(y)!=len(x):
        k[1] = f(x[len(y) - 1], y[len(y) - 1])
        k[2] = f(x[len(y) - 1] + 1/2*h, y[len(y) - 1] + h*k[1]/2)
        k[3] = f(x[len(y) - 1] + 1/2*h, y[len(y) - 1] + h*k[2]/2)
        k[4] = f(x[len(y) - 1] + h, y[len(y) - 1] + h*k[3])
        y.insert(len(y), y[len(y) - 1] +h*(k[1] + 2*k[2] + 2*k[3] + k[4])/6)
    return y




print ("Введите кол-во шагов")
h = int(input())
a = 0

print ("Введите конец интервала")
b = int(input())

print ("Введите точность")
e = float(input())
y1 = array('f')
x, step = np.linspace(a, b, h, True, True)
y1 = -x*x + 2*x - 2 + 12*math.e**(-x)
y = array('f', [10])
yA = array('f', [10])
yRKM = array('f', [10])


for i in x:
    if i != x[len(x)-1]:
        y.insert(len(y), Euler(y, x, step))

yA = Adams(x,yA,5,step)
xRKM,yRKM = RRKM(b, step, e)
yRK4 = RK4(x, step, b, h)



fig, ax = plt.subplots()
ax.set_xlim(0,5)
ax.set_ylim(-20,20)
ax.set_title("Графики")
ax.set_xlabel("x")
ax.set_ylabel("y")
ax.plot(x, y, color = 'g', label = 'Euler')
ax.plot(x, yRK4, color = 'y', label = 'RK4')
ax.plot(x, yA, color = 'c', label = 'Adams')
ax.plot(xRKM,yRKM, color = 'm', label = 'RKM')
ax.scatter(xRKM,yRKM,c = 'deeppink')
ax.plot(x, y1, color = 'r', label = 'Solution')
ax.legend()
plt.show()