import numpy as np
import math
import matplotlib.pyplot as plt

def f(x):
    return x - math.sin(x) - 0.25


a = -1
b = 1

n = 6
xr = np.linspace(a, b, n);
yr = np.array(list(map(f, xr)))

xch = np.zeros(n)

for i in range(n):
    xch[i] = 0.5*((b- a)*math.cos((2*i + 1) / (2*n ) * math.pi) + (b + a))

ych = np.array(list(map(f, xch)))

xx = xr


yy = yr

def lagrange(x):
    def l(k):
        a = 1
        for i in range(n):
            if i != k:
                a *= (x - xx[i])
        b = 1
        for i in range(n):
            if i != k:
                b *= (xx[k] - xx[i])
        return a / b;
    res = 0
    for i in range(n):
        res += yy[i]*l(i)
    return res





xg = np.arange(a, b + 0.1, 0.1)

maxr = 0

for cur in xg:
    maxr = max(maxr, abs(f(cur) - lagrange(cur)))

print(maxr)



def f1(x):
    return 1 - math.cos(x)

def f2(x):
    return math.sin(x)

p = 1
k = 3

xe = np.linspace(a, b, k)


ye0 = np.array(list(map(f, xe)))
ye1 = np.array(list(map(f1, xe)))
ye2 = np.array(list(map(f2, xe)))

m = k*(p + 1)

X = np.zeros((m, m))
for i in range(k):
    for j in range(p + 1):
        for l in range(0, m - j):
            X[i*(p+ 1) + j][l] = (xe[i]**(m -  1 - l - j)) * math.factorial(m - 1- l)/math.factorial(m - 1- j - l)

Y = np.zeros(m)

Y[0] = ye0[0]
Y[1] = ye1[0]

Y[2] = ye0[1]
Y[3] = ye1[1]

Y[4] = ye0[2]
Y[5] = ye1[2]


h = np.linalg.solve(X, Y)

print(X)
print(Y)

def e(x):
    res = 0
    for i in range(m):
        res += (x**i) * (h[m - 1- i])
    return res

def dedx(x):
    res = 0
    for i in range(m):
        res += i*(x**(i - 1)) * (h[m - 1- i])
    return res

xg = np.arange(a, b + 0.1, 0.1)

maxr = 0;
for curx in xg:
    maxr = max(maxr, abs(f(curx) - e(curx)))


yg = np.array(list(map(f1, xg)))
plt.plot(xg, yg)


yg = np.array(list(map(dedx, xg)))
plt.plot(xg, yg)

p = 0
k = 6

xe = np.linspace(a, b, k)

m = k*(p + 1)

X = np.zeros((m, m))
for i in range(k):
    for j in range(p + 1):
        for l in range(0, m - j):
            X[i*(p+ 1) + j][l] = (xe[i]**(m -  1 - l - j)) * math.factorial(m - 1- l)/math.factorial(m - 1- j - l)

Y = np.zeros(m)

ye0 = np.array(list(map(f, xe)))

Y[0] = ye0[0]
Y[1] = ye0[1]

Y[2] = ye0[2]
Y[3] = ye0[3]

Y[4] = ye0[4]
Y[5] = ye0[5]


h = np.linalg.solve(X, Y)



yg = np.array(list(map(dedx, xg)))
plt.plot(xg, yg)




plt.show()