matbiz01

Untitled

Apr 29th, 2021
610
Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
  1. import matplotlib.pyplot as plt
  2.  
  3. from itertools import zip_longest
  4. from math import sin, cos
  5. from math import pi
  6.  
  7. #tutaj możemy ustawić ilość węzłów - n, oraz sposób dobrania węzłów (czebyszew = 0 -> węzły równoodległe)
  8. x0 = -1
  9. xn = 1
  10. n = 10
  11. czebyszew = 0
  12.  
  13. def getOrigValues(x):
  14.     return 1 / ( 1 + 25 * x * x)
  15.  
  16.  
  17. def calculateBasicPoints(n, czebyszew):
  18.     points = []
  19.     interval = (xn - x0) / (n - 1)
  20.     for x in range(n):
  21.         if czebyszew:
  22.             currX = 1/2 * (xn + x0) + 1/2 * (xn - x0) * cos((2 * (x + 1) - 1) / 2 / n * pi)
  23.         else:
  24.             currX = x0 + x * interval
  25.         points.append([currX, getOrigValues(currX)])
  26.     return points
  27.  
  28.  
  29. #kod dla metody lagrange'a z użyciem wzoru newtona
  30. #aby go uruchomić, należy go odkomentować, i zakomentować kod dla wzoru newtona
  31.  
  32.  
  33. def makePyramid(basePoints, n):
  34.     xVals = [x for x,y in basePoints]
  35.     pyramid = []
  36.     pyramid.append(basePoints[0])
  37.     pyramid.append(basePoints[1])
  38.     pyramid[1].append((basePoints[1][1] - basePoints[0][1]) / (basePoints[1][0] - basePoints[0][0]))
  39.     for x in range(2, n):
  40.         pyramid.append(basePoints[x])
  41.         for y in range(x):
  42.             pyramid[x].append((pyramid[x][y + 1] - pyramid[x - 1][y + 1]) / (xVals[x] - xVals[x - (y + 1)])) #no idea how, but this works
  43.     coefficients = []
  44.     for x in range(n):
  45.         coefficients.append(pyramid[x][x + 1])
  46.     return coefficients
  47.  
  48.  
  49. def getNewtonVal(coefficients, x, xVals):
  50.     points = list(reversed(xVals))
  51.     coeffs = list(reversed(coefficients))
  52.     total = coeffs[0]
  53.     for i in range(1, len(coeffs)):
  54.         total *= (x - points[i])
  55.         total += coeffs[i]
  56.     return total
  57.  
  58. def maxDev(l1, l2):
  59.     temp = 0
  60.     for x,y in zip(l1, l2):
  61.         if abs(x - y) > temp: temp = abs(x - y)
  62.     return temp
  63.  
  64.  
  65. def do(n):
  66.     basePoints = calculateBasicPoints(n, czebyszew)
  67.     xVals = [x for x,y in basePoints]
  68.     coefficients = makePyramid(basePoints, n)
  69.  
  70.     point_count = 500
  71.     pointX = []
  72.     pointY = []
  73.     pointYorig = []
  74.     interval = (xn - x0) / point_count
  75.     for x in range(point_count):
  76.         currX = x0 + x * interval
  77.         pointX.append(currX)
  78.         pointY.append(getNewtonVal(coefficients, currX, xVals))
  79.         pointYorig.append(getOrigValues(currX))
  80.  
  81.     print(n, maxDev(pointY, pointYorig))
  82.     plt.plot(pointX, pointY, label="wielomian interpolujący")
  83.     plt.plot(pointX, pointYorig, label="oryginalna funkcja")
  84.     for point in calculateBasicPoints(n, czebyszew):
  85.         plt.scatter(point[0], point[1], s=50)
  86.     plt.legend(loc=1)
  87.     plt.show()
  88.    
  89.  
  90.  
  91. do(n)
  92.  
  93.  
  94.  
  95. #kod dla metody lagrange'a z użyciem metody lagrange'a
  96.  
  97. '''
  98. class LarangeTerm:
  99.    def __init__(self, yCoeff, denominator, xVals):
  100.        self.yCoeff  = yCoeff
  101.        self.denominator = denominator
  102.        self.xVals = xVals
  103.    def getValAt(self, x):
  104.        total = self.yCoeff / self.denominator
  105.        for xVal in self.xVals:
  106.            total *= (x - xVal)
  107.        return total
  108.    def print(self):
  109.        toPrint = ""
  110.        toPrint += str(self.yCoeff)
  111.        for x in self.xVals:
  112.            toPrint += "(x - " + str(x) + ")"
  113.        toPrint += "/" + str(self.denominator)
  114.        print(toPrint)
  115.  
  116. def lagrangeDenominators(z, xVals):
  117.    denominators = []
  118.    for n in range(z):
  119.        baseX = xVals[n]
  120.        total = 1
  121.        for i in range(len(xVals)):
  122.            if n != i:
  123.                total *= (baseX - xVals[i])
  124.        denominators.append(total)
  125.    return denominators
  126.  
  127. def makeLarTerms(n):
  128.    xVals = [x for x, y in calculateBasicPoints(n, czebyszew)]
  129.    yVals = [y for x, y in calculateBasicPoints(n, czebyszew)]
  130.    LarTerms = []
  131.    lagrDenoms = lagrangeDenominators(n, xVals)
  132.    xToLar = []
  133.    for k in range(n):
  134.        for i in range(n):
  135.            if i != k:
  136.                xToLar.append(xVals[i])
  137.        LarTerms.append(LarangeTerm(yVals[k], lagrDenoms[k], xToLar))
  138.        xToLar = []
  139.    return LarTerms
  140.  
  141. def lagrValAt(x, larTerms):
  142.    total = 0
  143.    for term in larTerms:
  144.        total += term.getValAt(x)
  145.    return total
  146.  
  147.  
  148. def maxDev(l1, l2):
  149.    temp = 0
  150.    for x,y in zip(l1, l2):
  151.        if abs(x - y) > temp: temp = abs(x - y)
  152.    return temp
  153.  
  154. def draw(point_count, n):
  155.    lart = makeLarTerms(n)
  156.  
  157.    pointX = []
  158.    pointY = []
  159.    pointYorig = []
  160.    interval = (xn - x0) / point_count
  161.    for x in range(point_count):
  162.        currX = x0 + x * interval
  163.        pointX.append(currX)
  164.        pointY.append(lagrValAt(currX, lart))
  165.        pointYorig.append(getOrigValues(currX))
  166.  
  167.  
  168.  
  169.  
  170.    print(n ,maxDev(pointY, pointYorig))
  171.  
  172.    
  173.    plt.plot(pointX, pointY, label="wielomian interpolujący")
  174.    plt.plot(pointX, pointYorig, label="oryginalna funkcja")
  175.    for point in calculateBasicPoints(n, czebyszew):
  176.        plt.scatter(point[0], point[1], s=50)
  177.    plt.legend(loc=1)
  178.    plt.show()
  179.  
  180.    
  181.  
  182. point_count = 1000
  183. draw(point_count, n)
  184.  
  185.  
  186. '''
  187.  
  188.  
  189.  
RAW Paste Data

Adblocker detected! Please consider disabling it...

We've detected AdBlock Plus or some other adblocking software preventing Pastebin.com from fully loading.

We don't have any obnoxious sound, or popup ads, we actively block these annoying types of ads!

Please add Pastebin.com to your ad blocker whitelist or disable your adblocking software.

×