ENH: optimize get_value_opt in class Function

rocketpy/mathutils/function.py

@@ -488,19 +488,18 @@ def get_value_opt(x):
 
             # change the function's name to avoid mypy's error
             def get_value_opt_multiple(*args):
-                x = np.array([[float(x) for x in list(args)]])
-                numerator_sum = 0
-                denominator_sum = 0
-                for i in range(len_y_data):
-                    sub = x_data[i] - x
-                    distance = np.linalg.norm(sub)
-                    if distance == 0:
-                        numerator_sum = y_data[i]
-                        denominator_sum = 1
-                        break
-                    weight = distance ** (-3)
-                    numerator_sum = numerator_sum + y_data[i] * weight
-                    denominator_sum = denominator_sum + weight
+                x = np.array([[float(val) for val in args]])
+                sub_matrix = x_data - x
+                distances_squared = np.sum(sub_matrix**2, axis=1)
+
+                zero_distance_index = np.where(distances_squared == 0)[0]
+                if len(zero_distance_index) > 0:
+                    return y_data[zero_distance_index[0]]
+
+                weights = distances_squared ** (-1.5)
+                numerator_sum = np.sum(y_data * weights)
+                denominator_sum = np.sum(weights)
+
                 return numerator_sum / denominator_sum
 
             get_value_opt = get_value_opt_multiple

tests/unit/test_function.py

@@ -123,3 +123,38 @@ def test_integral_function():
     """
     zero_func = Function(0)
     assert isinstance(zero_func, Function)
+
+
+@pytest.mark.parametrize(
+    "x, y, z",
+    [
+        (0, 0, 11.22540929),
+        (1, 2, 10),
+        (2, 3, 15.272727273),
+        (3, 4, 20),
+        (4, 5, 24.727272727),
+        (5, 6, 30),
+        (10, 10, 25.7201184),
+    ],
+)
+def test_get_value_opt(x, y, z):
+    """Test the get_value_opt method of the Function class. Currently only tests
+    a 2D function with shepard interpolation method. The tests are made on
+    points that are present or not in the source array.
+
+    Parameters
+    ----------
+    x : scalar
+        The x-coordinate.
+    y : scalar
+        The y-coordinate.
+    z : float
+        The expected interpolated value at (x, y).
+    """
+    x_data = np.array([1.0, 3.0, 5.0])
+    y_data = np.array([2.0, 4.0, 6.0])
+    z_data = np.array([10.0, 20.0, 30.0])
+    source = np.column_stack((x_data, y_data, z_data))
+    func = Function(source, interpolation="shepard", extrapolation="natural")
+    assert isinstance(func.get_value_opt(x, y), float)
+    assert np.isclose(func.get_value_opt(x, y), z, atol=1e-6)