```python
# -*- coding: utf-8 -*-
"""
Created on Thu Sep 16 20:23:07 2021

@author: gansa001
"""
from sympy import *
from sympy.plotting import plot
import matplotlib.pyplot as plt
import numpy as np
import math
```

### 1. Write a computer program to calculate the Lagrange interpolation polynomial Pn(x) to f(x) such that: Pn(xj)=f(xj)  j=0,1,2,……,n



```python
def delta(j,a,x,N):
    '''
    Parameters
    ----------
    j : int. 
        the jth lagrange polynomial
    a : sympy symbol
        the dependent term in the lagrange polynomial
    x : list
        x values from the original function
    N : int

    Returns
    -------
    The lagrange polynomial for any given j

    '''
    num=1
    den=1
    for i in range(N):
        if i!=j:
            num*=(a-x[i])
            den*=(x[j]-x[i])
    return (num/den)

def Lagrange(x,y):
    '''

    Parameters
    ----------
    x : list
        dependent variable values
    y : list
        independent variable values

    Returns
    -------
    f: The lagrangian approximation of the function

    '''
    N=len(x)
    a=Symbol('a')
    f=0
    for j in range(len(x)):
        answer= (delta(j,a,x,N))
        f+=y[j]*answer
    return (simplify(f))
```

#### a)	Apply it to some elementary functions or discrete data sets.


```python
'''
Exponential function from range -1 to 1 and 3 points to confirm the values from class. 
'''
x=[-1,0,1]
y=[math.exp(i) for i in x]

print('x = ',x)
print('y = ',y)
```

    x =  [-1, 0, 1]
    y =  [0.36787944117144233, 1.0, 2.718281828459045]



```python
##function call
Lagrange(x,y)
```




$\displaystyle 0.543080634815244 a^{2} + 1.1752011936438 a + 1.0$




```python
##plotting
a=Symbol('a')

p=plot(Lagrange(x,y), (E**a),(a,-1,1),show=False)
p[1].line_color='r'
p[0].label='Pn(x)'
p[1].label='f(x)'
p.legend=True
p.show()
```


```python
'''
Quartic function from range -1 to 1 and 3 points. 
'''
x=[-1,0,1]
y=[i**4 for i in x]

print('x = ',x)
print('y = ',y)
```

    x =  [-1, 0, 1]
    y =  [1, 0, 1]



```python
Lagrange(x,y)
```




$\displaystyle a^{2}$




```python
p=plot(Lagrange(x,y), (a**4),(a,-1,1),show=False)
p[1].line_color='r'
p[0].label='Pn(x)'
p[1].label='f(x)'
p.legend=True
p.show()
```

#### b)	Apply it to function $$ f(x) = \frac{1}{1+25x^2} $$  in the interval [-1,1] with equally spaced interpolation points and n=8, 16, 32.  Plot Pn(x) and f(x) in [-1,1] interval to see the interpolation error. 


```python
n=[8,16,32]

for npoints in n:
    print('********** NUMBER OF POINTS = ',npoints,' **********')
    print()
    x=list(np.linspace(-1,1,npoints))
    y=[(1/(1+25*(i**2))) for i in x]
    
    print('Lagrange= ',Lagrange(x,y))
    
    p=plot(Lagrange(x,y), (1/(1+25*(a**2))),(a,-1,1),show=False)
    p[1].line_color='r'
    p[0].label='Pn(x)'
    p[1].label='f(x)'
    p.legend=True
    p.show()
    
    print('\n\n')
    
```

#### c)	Do the same as in b, but use the roots of the Chebychev polynomial Tn+1(x) as interpolation points. Compare results with those obtained in (b).


```python
for npoints in n:
    print('********** NUMBER OF POINTS = ',npoints,' **********')
    x = [math.cos(((math.pi/2)*(2*j+1))/npoints) for j in range(npoints)]
    print('Chebyshev roots= ',x)
    print()
    y=[(1/(1+25*(i**2))) for i in x]
    
    print('Lagrange= ',Lagrange(x,y))
    
    p=plot(Lagrange(x,y), (1/(1+25*(a**2))),(a,-1,1),show=False)
    p[1].line_color='r'
    p[0].label='Pn(x)'
    p[1].label='f(x)'
    p.legend=True
    p.show()
    
    print('\n\n')
```

### Comparing errors


```python
for npoints in n:
    print('********** NUMBER OF POINTS = ',npoints,' **********')
    x = [math.cos(((math.pi/2)*(2*j+1))/npoints) for j in range(npoints)]
    print('Chebyshev roots= ',x)
    print()
    x1=list(np.linspace(-1,1,npoints))
    y=[(1/(1+25*(i**2))) for i in x]
    
    print('Lagrange= ',Lagrange(x,y))
    
    p=plot(Lagrange(x,y),Lagrange(x1,y), (1/(1+25*(a**2))),(a,-1,1),show=False)
    p[1].line_color='r'
    p[2].line_color='m'
    p[1].label='Pn(x)'
    p[0].label='Tn(x)'
    p[2].label='f(x)'
    p.legend=True
    p.show()
    
    print('\n\n')
```

# The plots above show that using the Chebyshev roots provide a much smaller error when doing Lagragian interpolation as opposed to using equally spaced points in a given range.

0cc319a6271dbe222c3b607104b7c90822e65a1b
ipynb
Jupyter Notebook
Lagrange-Chebyshev Interpolation Error.ipynb
GJAnsah/Lagrangian
8619f905fff0943242e3069404f49d55d8ca3f5a
MIT
Lagrange-Chebyshev Interpolation Error.ipynb
GJAnsah/Lagrangian
8619f905fff0943242e3069404f49d55d8ca3f5a
MIT
Lagrange-Chebyshev Interpolation Error.ipynb
GJAnsah/Lagrangian
8619f905fff0943242e3069404f49d55d8ca3f5a
MIT
Qwen/Qwen-72B
1. YES
2. YES