matbiz01

newton

Apr 8th, 2021
443
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.  
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.

×