Interpolation Exercise 2


In [1]:
%matplotlib inline
import matplotlib.pyplot as plt
import seaborn as sns
import numpy as np
sns.set_style('white')

In [2]:
from scipy.interpolate import griddata

Sparse 2d interpolation

In this example the values of a scalar field $f(x,y)$ are known at a very limited set of points in a square domain:

  • The square domain covers the region $x\in[-5,5]$ and $y\in[-5,5]$.
  • The values of $f(x,y)$ are zero on the boundary of the square at integer spaced points.
  • The value of $f$ is known at a single interior point: $f(0,0)=1.0$.
  • The function $f$ is not known at any other points.

Create arrays x, y, f:

  • x should be a 1d array of the x coordinates on the boundary and the 1 interior point.
  • y should be a 1d array of the y coordinates on the boundary and the 1 interior point.
  • f should be a 1d array of the values of f at the corresponding x and y coordinates.

You might find that np.hstack is helpful.


In [3]:
#xx= np.array([[0],[-5],[-4],[-3],[-2],[-1],[0],[1],[2],[3],[4],[5],[-5],[-4],[-3],[-2],[-1],[0],[1],[2],[3],[4],[5],[-5],[-5],[-5],[-5],[-5],[-5],[-5],[-5],[-5],[5],[5],[5],[5],[5],[5],[5],[5],[5]])
#yy= np.array([[0],[-5],[-5],[-5],[-5],[-5],[-5],[-5],[-5],[-5],[-5],[-5],[5],[5],[5],[5],[5],[5],[5],[5],[5],[5],[5],[-4],[-3],[-2],[-1],[0],[1],[2],[3],[4],[-4],[-3],[-2],[-1],[0],[1],[2],[3],[4]])
x= np.array([0,-5,-4,-3,-2,-1,0,1,2,3,4,5,-5,-4,-3,-2,-1,0,1,2,3,4,5,-5,-5,-5,-5,-5,-5,-5,-5,-5,5,5,5,5,5,5,5,5,5])
y=np.array([0,-5,-5,-5,-5,-5,-5,-5,-5,-5,-5,-5,5,5,5,5,5,5,5,5,5,5,5,-4,-3,-2,-1,0,1,2,3,4,-4,-3,-2,-1,0,1,2,3,4])
#f=np.array([1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0])

def F(x,y):
    f=[]
    for i in range(len(x)):               #Noah Miller helped me out with this part
        if abs(x[i])==5 or abs(y[i])==5:
            f.append(0)
        elif x[i]==0 and y[i]==0:
            f.append(1)
    return np.array(f)

f=F(x,y)
f


Out[3]:
array([1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
       0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0])

The following plot should show the points on the boundary and the single point in the interior:


In [4]:
plt.scatter(x, y);



In [5]:
assert x.shape==(41,)
assert y.shape==(41,)
assert f.shape==(41,)
assert np.count_nonzero(f)==1

Use meshgrid and griddata to interpolate the function $f(x,y)$ on the entire square domain:

  • xnew and ynew should be 1d arrays with 100 points between $[-5,5]$.
  • Xnew and Ynew should be 2d versions of xnew and ynew created by meshgrid.
  • Fnew should be a 2d array with the interpolated values of $f(x,y)$ at the points (Xnew,Ynew).
  • Use cubic spline interpolation.

In [10]:
"""points : ndarray of floats, shape (n, D)
    Data point coordinates. Can either be an array of
    shape (n, D), or a tuple of `ndim` arrays.
values : ndarray of float or complex, shape (n,)
    Data values.
xi : ndarray of float, shape (M, D)
    Points at which to interpolate data."""
#np.meshgrid?
#griddata(points,values,xi,method='linear')

xnew=np.linspace(-5,5,100)
ynew=np.linspace(-5,5,100)

Xnew, Ynew = np.meshgrid(xnew,ynew)
Fnew=F(Xnew,Ynew)


---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
<ipython-input-10-0dbf6b406fae> in <module>()
     13 
     14 Xnew, Ynew = np.meshgrid(xnew,ynew)
---> 15 Fnew=F(Xnew,Ynew)

<ipython-input-3-f74075e6eda8> in F(x, y)
      9     f=[]
     10     for i in range(len(x)):               #Noah Miller helped me out with this part
---> 11         if abs(x[i])==5 or abs(y[i])==5:
     12             f.append(0)
     13         elif x[i]==0 and y[i]==0:

ValueError: The truth value of an array with more than one element is ambiguous. Use a.any() or a.all()

In [11]:
assert xnew.shape==(100,)
assert ynew.shape==(100,)
assert Xnew.shape==(100,100)
assert Ynew.shape==(100,100)
assert Fnew.shape==(100,100)

Plot the values of the interpolated scalar field using a contour plot. Customize your plot to make it effective and beautiful.


In [12]:
plt.pcolor(Xnew, Ynew, Fnew);
plt.colorbar();
plt.scatter(Xnew, Ynew, marker='o', color='blue', label='interpolated points')
plt.title('Interpolated Scalar Field')
plt.xlabel('x')
plt.ylabel('y');



In [ ]:
assert True # leave this to grade the plot