In [2]:
%matplotlib inline
import matplotlib.pyplot as plt
import seaborn as sns
import numpy as np
sns.set_style('white')
In [3]:
from scipy.interpolate import griddata
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:
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 [6]:
x1=np.arange(-5,6)
y1=5*np.ones(11)
f1=np.zeros(11)
x2=np.arange(-5,6)
y2=-5*np.ones(11)
f2=np.zeros(11)
y3=np.arange(-4,5)
x3=5*np.ones(9)
f3=np.zeros(9)
y4=np.arange(-4,5)
x4=-5*np.ones(9)
f4=np.zeros(9)
x5=np.array([0])
y5=np.array([0])
f5=np.array([1])
x=np.hstack((x1,x2,x3,x4,x5))
y=np.hstack((y1,y2,y3,y4,y5))
f=np.hstack((f1,f2,f3,f4,f5))
print (x)
print (y)
print (f)
The following plot should show the points on the boundary and the single point in the interior:
In [7]:
plt.scatter(x, y);
In [8]:
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).
In [10]:
xnew=np.linspace(-5.0,6.0,100)
ynew=np.linspace(-5,6,100)
Xnew,Ynew=np.meshgrid(xnew,ynew)
Fnew=griddata((x,y),f,(Xnew,Ynew),method='cubic',fill_value=0.0)
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 [16]:
plt.figure(figsize=(10,8))
plt.contourf(Fnew,cmap='cubehelix_r')
plt.title('2D Interpolation');
In [ ]:
assert True # leave this to grade the plot