added econforge interpolation into LinearInterp and scipy interpolati…

HARK/fast_interpolation.py

@@ -0,0 +1,16 @@
+import numpy as np
+from interpolation import interp # econ forge
+from scipy.interpolate import CubicSpline # scipy
+
+# econforge interp and then uses the hark extrapolation for values outside of the range
+def fast_linear_1d_interp(x, obj):
+    result_temp = interp(obj.x_list, obj.y_list, x)
+    above_upper_bound = x > obj.x_list[-1]
+    i = len(obj.x_list) - 1
+    alpha = (x[above_upper_bound] - obj.x_list[i - 1]) / (obj.x_list[i] - obj.x_list[i - 1])
+    result_temp[above_upper_bound] = (1.0 - alpha) * obj.y_list[i - 1] + alpha * obj.y_list[i]
+    return [result_temp]
+
+# gets the scipy cubic spline object
+def get_cubic_spline(obj):
+    return CubicSpline(np.array(obj.x_list), np.array(obj.y_list))

HARK/interpolation.py

@@ -15,6 +15,9 @@
 from HARK.utilities import CRRAutility, CRRAutilityP, CRRAutilityPP
 from .core import MetricObject
 
+from HARK.fast_interpolation import fast_linear_1d_interp, get_cubic_spline
+from interpolation.splines import CGrid, eval_linear
+
 
 def _isscalar(x):
     """
@@ -814,6 +817,15 @@ def _evalOrDer(self, x, _eval, _Der):
         A list including the level and/or derivative of the interpolated function where requested.
         """
 
+        # ========== fast interpolation check for econ hark library ==========
+        # sometimes x is empty causing the below to error, so added a try except
+        # tested with an extra check to see if x is empty, but this added to the run time
+        try:
+            if _eval and not _Der and not self.decay_extrap and x[0] >= self.x_list[0]:
+                return fast_linear_1d_interp(x, self)
+        except:
+            print("empty x sent to LinearInterp")
+
         i = np.maximum(np.searchsorted(self.x_list[:-1], x), 1)
         alpha = (x - self.x_list[i - 1]) / (self.x_list[i] - self.x_list[i - 1])
 
@@ -1014,6 +1026,7 @@ def _evaluate(self, x):
                     - self.coeffs[pos, 2] * np.exp(alpha * self.coeffs[pos, 3])
                 )
         else:
+
             m = len(x)
             pos = np.searchsorted(self.x_list, x)
             y = np.zeros(m)
@@ -1374,7 +1387,15 @@ def __init__(self, f_values, x_list, y_list, xSearchFunc=None, ySearchFunc=None)
         self.xSearchFunc = xSearchFunc
         self.ySearchFunc = ySearchFunc
 
+        self.customgrid = CGrid(self.x_list, self.y_list)
+
     def _evaluate(self, x, y):
+
+        try:
+            return eval_linear(self.customgrid, self.f_values, np.column_stack((x,y)))
+        except:
+            print("econforge point format not supported : ",np.column_stack((x,y)))
+
         """
         Returns the level of the interpolated function at each value in x,y.
         Only called internally by HARKinterpolator2D.__call__ (etc).
@@ -1532,7 +1553,15 @@ def __init__(
         self.ySearchFunc = ySearchFunc
         self.zSearchFunc = zSearchFunc
 
+        self.customgrid = CGrid(self.x_list, self.y_list, self.z_list)
+
     def _evaluate(self, x, y, z):
+
+        try:
+            return eval_linear(self.customgrid, self.f_values, np.column_stack((x,y,z)))
+        except:
+            print("econforge point format not supported : ",np.column_stack((x,y,z)))
+
         """
         Returns the level of the interpolated function at each value in x,y,z.
         Only called internally by HARKinterpolator3D.__call__ (etc).