modify akima to not allocate on repeated calls

byuflowlab · Jan 24, 2020 · eedf106 · eedf106
1 parent ad81cea
commit eedf106
Show file tree

Hide file tree

Showing 3 changed files with 14 additions and 8 deletions.
diff --git a/README.md b/README.md
@@ -65,8 +65,9 @@ or if you plan to evaluate the spline repeatedly
 ```julia
 
 spline = akima_setup(x, y)
-ypt1 = akima_interp(xpt1, spline)
-ypt2 = akima_interp(xpt2, spline)
+ypt = similar(xpt1)
+akima_interp!(ypt, xpt1, spline)  # ypt changed in place
+akima_interp!(ypt, xpt2, spline)
 ```
 
 ### 2D/3D/4D Interpolation using Recursive 1D Interpolation

diff --git a/src/interpolate.jl b/src/interpolate.jl
@@ -87,16 +87,15 @@ end
 
 
 """
-    akima_interp(x, spline)
+    akima_interp!(y, x, spline)
 
 Evaluate the spline setup by akima_setup.  x are the locations to evaluate the spline at
-and spline is the object returned by akima_setup.
+and spline is the object returned by akima_setup.  the result is returned in place at y.
 """
-function akima_interp(x, spline)
+function akima_interp!(y, x, spline)
 
     # interpolate at each point
     n = length(x)
-    y = zeros(eltype(x), n)
 
     for i = 1:n
 
@@ -111,7 +110,7 @@ function akima_interp(x, spline)
 
     end
 
-    return y
+    return nothing
 end
 
 """
@@ -122,7 +121,9 @@ xpt may be an array.
 """
 function akima(xdata, ydata, xpt)
     spline = akima_setup(xdata, ydata)
-    return akima_interp(xpt, spline)
+    ypt = similar(xpt)
+    akima_interp!(ypt, xpt, spline)
+    return ypt
 end
 
 

diff --git a/test/runtests.jl b/test/runtests.jl
@@ -490,6 +490,10 @@ xpt = 0:pi/16:2*pi
 
 ypt = akima(x, y, xpt)
 
+# FIXME: this is not an independent test.  right now it is just testing against refactoring.
+ytest = [0.0, 0.21394356492851746, 0.40088834764831843, 0.5641656518405971, 0.7071067811865475, 0.8281808021479407, 0.920495128834866, 0.9793385864009951, 1.0, 0.9873810649285174, 0.9419417382415921, 0.8523082377305078, 0.7071067811865476, 0.5303300858899105, 0.35355339059327384, 0.17677669529663714, 1.2246467991473532e-16, -0.17677669529663698, -0.35355339059327373, -0.5303300858899104, -0.7071067811865475, -0.8523082377305078, -0.941941738241592, -0.9873810649285175, -1.0, -0.9793385864009951, -0.9204951288348657, -0.8281808021479408, -0.7071067811865477, -0.5641656518405971, -0.400888347648319, -0.21394356492851788, 0.0]
+
+@test isapprox(ypt, ytest, atol=1e-8)
 # ---------------------------