import matplotlib.pyplot as plt

from itertools import zip_longest
from math import sin, cos
from math import pi

#tutaj możemy ustawić ilość węzłów - n, oraz sposób dobrania węzłów (czebyszew = 0 -> węzły równoodległe)
x0 = -1
xn = 1
n = 10
czebyszew = 0

def getOrigValues(x):
    return 1 / ( 1 + 25 * x * x)


def calculateBasicPoints(n, czebyszew):
    points = []
    interval = (xn - x0) / (n - 1)
    for x in range(n):
        if czebyszew:
            currX = 1/2 * (xn + x0) + 1/2 * (xn - x0) * cos((2 * (x + 1) - 1) / 2 / n * pi)
        else:
            currX = x0 + x * interval
        points.append([currX, getOrigValues(currX)])
    return points


#kod dla metody lagrange'a z użyciem wzoru newtona
#aby go uruchomić, należy go odkomentować, i zakomentować kod dla wzoru newtona


def makePyramid(basePoints, n):
    xVals = [x for x,y in basePoints]
    pyramid = []
    pyramid.append(basePoints[0])
    pyramid.append(basePoints[1])
    pyramid[1].append((basePoints[1][1] - basePoints[0][1]) / (basePoints[1][0] - basePoints[0][0]))
    for x in range(2, n):
        pyramid.append(basePoints[x])
        for y in range(x):
            pyramid[x].append((pyramid[x][y + 1] - pyramid[x - 1][y + 1]) / (xVals[x] - xVals[x - (y + 1)])) #no idea how, but this works
    coefficients = []
    for x in range(n):
        coefficients.append(pyramid[x][x + 1])
    return coefficients


def getNewtonVal(coefficients, x, xVals):
    points = list(reversed(xVals))
    coeffs = list(reversed(coefficients))
    total = coeffs[0]
    for i in range(1, len(coeffs)):
        total *= (x - points[i])
        total += coeffs[i]
    return total

def maxDev(l1, l2):
    temp = 0
    for x,y in zip(l1, l2):
        if abs(x - y) > temp: temp = abs(x - y)
    return temp


def do(n):
    basePoints = calculateBasicPoints(n, czebyszew)
    xVals = [x for x,y in basePoints]
    coefficients = makePyramid(basePoints, n)

    point_count = 500
    pointX = []
    pointY = []
    pointYorig = []
    interval = (xn - x0) / point_count
    for x in range(point_count):
        currX = x0 + x * interval
        pointX.append(currX)
        pointY.append(getNewtonVal(coefficients, currX, xVals))
        pointYorig.append(getOrigValues(currX))

    print(n, maxDev(pointY, pointYorig))
    plt.plot(pointX, pointY, label="wielomian interpolujący")
    plt.plot(pointX, pointYorig, label="oryginalna funkcja")
    for point in calculateBasicPoints(n, czebyszew):
        plt.scatter(point[0], point[1], s=50)
    plt.legend(loc=1)
    plt.show()



do(n)



#kod dla metody lagrange'a z użyciem metody lagrange'a

'''
class LarangeTerm:
    def __init__(self, yCoeff, denominator, xVals):
        self.yCoeff  = yCoeff
        self.denominator = denominator
        self.xVals = xVals
    def getValAt(self, x):
        total = self.yCoeff / self.denominator
        for xVal in self.xVals:
            total *= (x - xVal)
        return total
    def print(self):
        toPrint = ""
        toPrint += str(self.yCoeff)
        for x in self.xVals:
            toPrint += "(x - " + str(x) + ")"
        toPrint += "/" + str(self.denominator)
        print(toPrint)

def lagrangeDenominators(z, xVals):
    denominators = []
    for n in range(z):
        baseX = xVals[n]
        total = 1
        for i in range(len(xVals)):
            if n != i:
                total *= (baseX - xVals[i])
        denominators.append(total)
    return denominators

def makeLarTerms(n):
    xVals = [x for x, y in calculateBasicPoints(n, czebyszew)]
    yVals = [y for x, y in calculateBasicPoints(n, czebyszew)]
    LarTerms = []
    lagrDenoms = lagrangeDenominators(n, xVals)
    xToLar = []
    for k in range(n):
        for i in range(n):
            if i != k:
                xToLar.append(xVals[i])
        LarTerms.append(LarangeTerm(yVals[k], lagrDenoms[k], xToLar))
        xToLar = []
    return LarTerms

def lagrValAt(x, larTerms):
    total = 0
    for term in larTerms:
        total += term.getValAt(x)
    return total


def maxDev(l1, l2):
    temp = 0
    for x,y in zip(l1, l2):
        if abs(x - y) > temp: temp = abs(x - y)
    return temp

def draw(point_count, n):
    lart = makeLarTerms(n)

    pointX = []
    pointY = []
    pointYorig = []
    interval = (xn - x0) / point_count
    for x in range(point_count):
        currX = x0 + x * interval
        pointX.append(currX)
        pointY.append(lagrValAt(currX, lart))
        pointYorig.append(getOrigValues(currX))




    print(n ,maxDev(pointY, pointYorig))


    plt.plot(pointX, pointY, label="wielomian interpolujący")
    plt.plot(pointX, pointYorig, label="oryginalna funkcja")
    for point in calculateBasicPoints(n, czebyszew):
        plt.scatter(point[0], point[1], s=50)
    plt.legend(loc=1)
    plt.show()



point_count = 1000
draw(point_count, n)


'''