HARK icon indicating copy to clipboard operation
HARK copied to clipboard

Cubic Interpolation on lower bound

Open alanlujan91 opened this issue 3 years ago • 2 comments

This is related to issue #971.

An issue with CubicInterp is that it previously did not work as expected for the lower bound of the range where the interpolation is defined. As an example:

from HARK.interpolation import CubicInterp
import numpy as np
f = lambda x: x**2
df = lambda x: 2*x
x = np.linspace(-1,1,5)
y = f(x)
dy = df(x)
interp = CubicInterp(x,y,dy)
interp(x)

Out[1]: array([ nan, 0.25, 0.  , 0.25, 1.  ])

Clearly, the first element should be 1.. This is fixed when you set lower_extrap=True.

interp = CubicInterp(x,y,dy, lower_extrap=True)
interp(x)

Out[2]: array([1.  , 0.25, 0.  , 0.25, 1.  ])

However, x=-1 is not extrapolation, as it is one of the data points given, so it should return a value without setting lower_extrap=True.

A previous attempt to resolving this was the following line in `CubicInterp':

#972

...
y[x == self.x_list.min()] = self.y_list.min()
...

but this only works for a particular set of cases. A better solution is

#1008

...
y[x == self.x_list[0]] = self.y_list[0]
... 

interp = CubicInterp(x,y,dy)
interp(x)

Out[2]: array([1.  , 0.25, 0.  , 0.25, 1.  ])

which explicitly plugs in the lower bound. At some point, this correction was overwritten again, so HARK code is back to #972.

Evidently, #972 does not break any tests, as it creeped back in without notice. The way I realized this is back on HARK is because it created a conflict on #1011.

alanlujan91 avatar Aug 19 '21 17:08 alanlujan91