Interpolation Exercise 2


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

In [165]:
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 [166]:
a = np.linspace(-5.0,5.0,11)

In [167]:
b = np.array([-5.0,5.0])

In [168]:
c = np.array([-5.0,0,5.0])

In [169]:
d = np.linspace(-5.0,-5.0,11)

In [170]:
g = np.array([-4.0,-4.0,-3.0,-3.0,-2.0,-2.0,-1.0,-1.0])

In [171]:
h = np.array([0,0,0])

In [172]:
i = np.array([1.0,1.0,2.0,2.0,3.0,3.0,4.0,4.0])

In [173]:
e = np.linspace(5.0,5.0,11)

In [174]:
x = np.hstack((a,b,b,b,b,c,b,b,b,b,a))

In [175]:
y = np.hstack((d,g,h,i,e))

In [176]:
f = np.zeros((41))

In [177]:
f[21] = 1

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


In [178]:
plt.scatter(x, y);
plt.grid(True)



In [179]:
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 [180]:
xnew = np.linspace(-5.0,5.0,100)
ynew = np.linspace(-5.0,5.0,100)

In [181]:
Xnew, Ynew = np.meshgrid(xnew,ynew)

In [182]:
Fnew = griddata((x,y), f, (Xnew, Ynew), method='cubic')

In [183]:
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 [184]:
f = plt.figure(figsize=(7,5))
plt.contour(Xnew,Ynew,Fnew,cmap='hsv')
plt.title('Contour Plot of Scalar Field')
plt.xlabel('x')
plt.ylabel('y')
plt.xlim(-5,5)
plt.ylim(-5,5)
plt.colorbar();



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