From 5a7f31523b54dacbbbbbca04fd37b00dd8a86166 Mon Sep 17 00:00:00 2001 From: NimaSarajpoor Date: Tue, 10 Feb 2026 21:34:09 -0500 Subject: [PATCH 01/11] improve njit_sdp --- stumpy/sdp.py | 7 +++++-- 1 file changed, 5 insertions(+), 2 deletions(-) diff --git a/stumpy/sdp.py b/stumpy/sdp.py index 551e8f2f2..b46e3c6f2 100644 --- a/stumpy/sdp.py +++ b/stumpy/sdp.py @@ -23,11 +23,14 @@ def _njit_sliding_dot_product(Q, T): out : numpy.ndarray Sliding dot product between `Q` and `T`. """ - m = Q.shape[0] + m = len(Q) l = T.shape[0] - m + 1 out = np.empty(l) for i in range(l): - out[i] = np.dot(Q, T[i : i + m]) + result = 0.0 + for j in range(m): + result += Q[j] * T[i + j] + out[i] = result return out From dd98c83eb87f83daa680f017b3767c0d5e774637 Mon Sep 17 00:00:00 2001 From: NimaSarajpoor Date: Tue, 10 Feb 2026 21:47:36 -0500 Subject: [PATCH 02/11] simplify function --- stumpy/sdp.py | 7 +------ 1 file changed, 1 insertion(+), 6 deletions(-) diff --git a/stumpy/sdp.py b/stumpy/sdp.py index b46e3c6f2..bd225503a 100644 --- a/stumpy/sdp.py +++ b/stumpy/sdp.py @@ -52,12 +52,7 @@ def _convolve_sliding_dot_product(Q, T): output : numpy.ndarray Sliding dot product between `Q` and `T`. """ - n = T.shape[0] - m = Q.shape[0] - Qr = np.flipud(Q) # Reverse/flip Q - QT = convolve(Qr, T) - - return QT.real[m - 1 : n] + return convolve(np.flipud(Q), T, mode="valid") def _sliding_dot_product(Q, T): From b581b3817ef5e523b57e8808bff8e83461d9ce3a Mon Sep 17 00:00:00 2001 From: NimaSarajpoor Date: Tue, 10 Feb 2026 22:20:37 -0500 Subject: [PATCH 03/11] add new sdp func and its test --- stumpy/sdp.py | 20 ++++++++++++++++++++ tests/test_sdp.py | 7 +++++++ 2 files changed, 27 insertions(+) diff --git a/stumpy/sdp.py b/stumpy/sdp.py index bd225503a..2dba7d0a9 100644 --- a/stumpy/sdp.py +++ b/stumpy/sdp.py @@ -55,6 +55,26 @@ def _convolve_sliding_dot_product(Q, T): return convolve(np.flipud(Q), T, mode="valid") +def _oaconvolve_sliding_dot_product(Q, T): + """ + Use scipy's oaconvolve to calculate the sliding dot product. + + Parameters + ---------- + Q : numpy.ndarray + Query array or subsequence + + T : numpy.ndarray + Time series or sequence + + Returns + ------- + output : numpy.ndarray + Sliding dot product between `Q` and `T`. + """ + return oaconvolve(np.ascontiguousarray(Q[::-1]), T, mode="valid") + + def _sliding_dot_product(Q, T): """ Compute the sliding dot product between `Q` and `T` diff --git a/tests/test_sdp.py b/tests/test_sdp.py index 41e956274..c0a8f6e2d 100644 --- a/tests/test_sdp.py +++ b/tests/test_sdp.py @@ -29,6 +29,13 @@ def test_convolve_sliding_dot_product(Q, T): npt.assert_almost_equal(ref_mp, comp_mp) +@pytest.mark.parametrize("Q, T", test_data) +def test_oaconvolve_sliding_dot_product(Q, T): + ref_mp = naive.rolling_window_dot_product(Q, T) + comp_mp = sdp._convolve_sliding_dot_product(Q, T) + npt.assert_almost_equal(ref_mp, comp_mp) + + @pytest.mark.parametrize("Q, T", test_data) def test_sliding_dot_product(Q, T): ref_mp = naive.rolling_window_dot_product(Q, T) From e9c0e2448bd4895d32668c697a0f051696d406a3 Mon Sep 17 00:00:00 2001 From: NimaSarajpoor Date: Tue, 10 Feb 2026 22:31:52 -0500 Subject: [PATCH 04/11] add pocketfft sdp and test --- stumpy/sdp.py | 32 ++++++++++++++++++++++++++++++++ tests/test_sdp.py | 7 +++++++ 2 files changed, 39 insertions(+) diff --git a/stumpy/sdp.py b/stumpy/sdp.py index 2dba7d0a9..a69cf68fc 100644 --- a/stumpy/sdp.py +++ b/stumpy/sdp.py @@ -75,6 +75,38 @@ def _oaconvolve_sliding_dot_product(Q, T): return oaconvolve(np.ascontiguousarray(Q[::-1]), T, mode="valid") +def _pocketfft_sliding_dot_product(Q, T): + """ + Use scipy.fft._pocketfft to compute + the sliding dot product. + + Parameters + ---------- + Q : numpy.ndarray + Query array or subsequence + + T : numpy.ndarray + Time series or sequence + + Returns + ------- + output : numpy.ndarray + Sliding dot product between `Q` and `T`. + """ + n = len(T) + m = len(Q) + next_fast_n = next_fast_len(n, real=True) + + tmp = np.empty((2, next_fast_n)) + tmp[0, :m] = Q[::-1] + tmp[0, m:] = 0.0 + tmp[1, :n] = T + tmp[1, n:] = 0.0 + fft_2d = r2c(True, tmp, axis=-1) + + return c2r(False, np.multiply(fft_2d[0], fft_2d[1]), n=next_fast_n)[m - 1 : n] + + def _sliding_dot_product(Q, T): """ Compute the sliding dot product between `Q` and `T` diff --git a/tests/test_sdp.py b/tests/test_sdp.py index c0a8f6e2d..774a69d7e 100644 --- a/tests/test_sdp.py +++ b/tests/test_sdp.py @@ -36,6 +36,13 @@ def test_oaconvolve_sliding_dot_product(Q, T): npt.assert_almost_equal(ref_mp, comp_mp) +@pytest.mark.parametrize("Q, T", test_data) +def test_pocketfft_sliding_dot_product(Q, T): + ref_mp = naive.rolling_window_dot_product(Q, T) + comp_mp = sdp._pocketfft_sliding_dot_product(Q, T) + npt.assert_almost_equal(ref_mp, comp_mp) + + @pytest.mark.parametrize("Q, T", test_data) def test_sliding_dot_product(Q, T): ref_mp = naive.rolling_window_dot_product(Q, T) From 6c14ab5c823050c72725ee48fa5b78dec8504eae Mon Sep 17 00:00:00 2001 From: NimaSarajpoor Date: Tue, 10 Feb 2026 22:33:28 -0500 Subject: [PATCH 05/11] fixed missing imports --- stumpy/sdp.py | 4 +++- 1 file changed, 3 insertions(+), 1 deletion(-) diff --git a/stumpy/sdp.py b/stumpy/sdp.py index a69cf68fc..57c628b48 100644 --- a/stumpy/sdp.py +++ b/stumpy/sdp.py @@ -1,6 +1,8 @@ import numpy as np from numba import njit -from scipy.signal import convolve +from scipy.fft import next_fast_len +from scipy.fft._pocketfft.basic import c2r, r2c +from scipy.signal import convolve, oaconvolve from . import config From 7c1e52c9829483a3e818faa2417292195292f3c9 Mon Sep 17 00:00:00 2001 From: NimaSarajpoor Date: Tue, 10 Feb 2026 23:25:12 -0500 Subject: [PATCH 06/11] fixed coverage --- tests/test_sdp.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/tests/test_sdp.py b/tests/test_sdp.py index 774a69d7e..785cace68 100644 --- a/tests/test_sdp.py +++ b/tests/test_sdp.py @@ -32,7 +32,7 @@ def test_convolve_sliding_dot_product(Q, T): @pytest.mark.parametrize("Q, T", test_data) def test_oaconvolve_sliding_dot_product(Q, T): ref_mp = naive.rolling_window_dot_product(Q, T) - comp_mp = sdp._convolve_sliding_dot_product(Q, T) + comp_mp = sdp._oaconvolve_sliding_dot_product(Q, T) npt.assert_almost_equal(ref_mp, comp_mp) From 98692857b68262eb9bebf42f91174c5541e57d64 Mon Sep 17 00:00:00 2001 From: NimaSarajpoor Date: Wed, 11 Feb 2026 00:38:34 -0500 Subject: [PATCH 07/11] refactored test functions --- tests/test_sdp.py | 47 +++++++++++++++++++---------------------------- 1 file changed, 19 insertions(+), 28 deletions(-) diff --git a/tests/test_sdp.py b/tests/test_sdp.py index 785cace68..369647761 100644 --- a/tests/test_sdp.py +++ b/tests/test_sdp.py @@ -1,3 +1,6 @@ +import inspect +import warnings + import naive import numpy as np import pytest @@ -15,36 +18,24 @@ ] -@pytest.mark.parametrize("Q, T", test_data) -def test_njit_sliding_dot_product(Q, T): - ref_mp = naive.rolling_window_dot_product(Q, T) - comp_mp = sdp._njit_sliding_dot_product(Q, T) - npt.assert_almost_equal(ref_mp, comp_mp) - - -@pytest.mark.parametrize("Q, T", test_data) -def test_convolve_sliding_dot_product(Q, T): - ref_mp = naive.rolling_window_dot_product(Q, T) - comp_mp = sdp._convolve_sliding_dot_product(Q, T) - npt.assert_almost_equal(ref_mp, comp_mp) - - -@pytest.mark.parametrize("Q, T", test_data) -def test_oaconvolve_sliding_dot_product(Q, T): - ref_mp = naive.rolling_window_dot_product(Q, T) - comp_mp = sdp._oaconvolve_sliding_dot_product(Q, T) - npt.assert_almost_equal(ref_mp, comp_mp) +def get_sdp_function_names(): + out = [] + for func_name, func in inspect.getmembers(sdp, inspect.isfunction): + if func_name.endswith("sliding_dot_product"): + out.append(func_name) - -@pytest.mark.parametrize("Q, T", test_data) -def test_pocketfft_sliding_dot_product(Q, T): - ref_mp = naive.rolling_window_dot_product(Q, T) - comp_mp = sdp._pocketfft_sliding_dot_product(Q, T) - npt.assert_almost_equal(ref_mp, comp_mp) + return out @pytest.mark.parametrize("Q, T", test_data) def test_sliding_dot_product(Q, T): - ref_mp = naive.rolling_window_dot_product(Q, T) - comp_mp = sdp._sliding_dot_product(Q, T) - npt.assert_almost_equal(ref_mp, comp_mp) + for func_name in get_sdp_function_names(): + func = getattr(sdp, func_name) + try: + comp = func(Q, T) + ref = naive.rolling_window_dot_product(Q, T) + npt.assert_allclose(comp, ref) + except Exception as e: # pragma: no cover + msg = f"Error in {func_name}, with n_Q={len(Q)} and n_T={len(T)}" + warnings.warn(msg) + raise e From 05b073c0b6358f47876055bf1f23e200f605cf35 Mon Sep 17 00:00:00 2001 From: NimaSarajpoor Date: Wed, 11 Feb 2026 00:44:59 -0500 Subject: [PATCH 08/11] minor change --- stumpy/sdp.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/stumpy/sdp.py b/stumpy/sdp.py index 57c628b48..bd1d8fea7 100644 --- a/stumpy/sdp.py +++ b/stumpy/sdp.py @@ -26,7 +26,7 @@ def _njit_sliding_dot_product(Q, T): Sliding dot product between `Q` and `T`. """ m = len(Q) - l = T.shape[0] - m + 1 + l = len(T) - m + 1 out = np.empty(l) for i in range(l): result = 0.0 From 8b6191beff61678b38fe21181958301198c90052 Mon Sep 17 00:00:00 2001 From: NimaSarajpoor Date: Wed, 11 Feb 2026 01:34:11 -0500 Subject: [PATCH 09/11] enhanced test cases --- tests/test_sdp.py | 136 ++++++++++++++++++++++++++++++++++++++++++---- 1 file changed, 124 insertions(+), 12 deletions(-) diff --git a/tests/test_sdp.py b/tests/test_sdp.py index 369647761..3f9b18949 100644 --- a/tests/test_sdp.py +++ b/tests/test_sdp.py @@ -1,5 +1,6 @@ import inspect import warnings +from operator import eq, lt import naive import numpy as np @@ -8,13 +9,61 @@ from stumpy import sdp -test_data = [ - (np.array([-1, 1, 2], dtype=np.float64), np.array(range(5), dtype=np.float64)), +# README +# Real FFT algorithm performs more efficiently when the length +# of the input array `arr` is composed of small prime factors. +# The next_fast_len(arr, real=True) function from Scipy returns +# the same length if len(arr) is composed of a subset of +# prime numbers 2, 3, 5. Therefore, these radices are +# considered as the most efficient for the real FFT algorithm. + +# To ensure that the tests cover different cases, the following cases +# are considered: +# 1. len(T) is even, and len(T) == next_fast_len(len(T), real=True) +# 2. len(T) is odd, and len(T) == next_fast_len(len(T), real=True) +# 3. len(T) is even, and len(T) < next_fast_len(len(T), real=True) +# 4. len(T) is odd, and len(T) < next_fast_len(len(T), real=True) +# And 5. a special case of 1, where len(T) is power of 2. + +# Therefore: +# 1. len(T) is composed of 2 and a subset of {3, 5} +# 2. len(T) is composed of a subset of {3, 5} +# 3. len(T) is composed of a subset of {7, 11, 13, ...} and 2 +# 4. len(T) is composed of a subset of {7, 11, 13, ...} +# 5. len(T) is power of 2 + +# In some cases, the prime factors are raised to a power of +# certain degree to increase the length of array to be around +# 1000-2000. This allows us to test sliding_dot_product for +# wider range of query lengths. + +test_inputs = [ + # Input format: + # ( + # len(T), + # remainder, # from `len(T) % 2` + # comparator, # for len(T) comparator next_fast_len(len(T), real=True) + # ) ( - np.array([9, 8100, -60], dtype=np.float64), - np.array([584, -11, 23, 79, 1001], dtype=np.float64), - ), - (np.random.uniform(-1000, 1000, [8]), np.random.uniform(-1000, 1000, [64])), + 2 * (3**2) * (5**3), + 0, + eq, + ), # = 2250, Even `len(T)`, and `len(T) == next_fast_len(len(T), real=True)` + ( + (3**2) * (5**3), + 1, + eq, + ), # = 1125, Odd `len(T)`, and `len(T) == next_fast_len(len(T), real=True)`. + ( + 2 * 7 * 11 * 13, + 0, + lt, + ), # = 2002, Even `len(T)`, and `len(T) < next_fast_len(len(T), real=True)` + ( + 7 * 11 * 13, + 1, + lt, + ), # = 1001, Odd `len(T)`, and `len(T) < next_fast_len(len(T), real=True)` ] @@ -27,15 +76,78 @@ def get_sdp_function_names(): return out -@pytest.mark.parametrize("Q, T", test_data) -def test_sliding_dot_product(Q, T): +@pytest.mark.parametrize("n_T, remainder, comparator", test_inputs) +def test_sdp(n_T, remainder, comparator): + # test_sdp for cases 1-4 + + n_Q_prime = [ + 2, + 3, + 5, + 7, + 11, + 13, + 17, + 19, + 23, + 29, + 31, + 37, + 41, + 43, + 47, + 53, + 59, + 61, + 67, + 71, + 73, + 79, + 83, + 89, + 97, + ] + n_Q_power2 = [2, 4, 8, 16, 32, 64] + n_Q_values = n_Q_prime + n_Q_power2 + [n_T] + n_Q_values = sorted(n_Q for n_Q in set(n_Q_values) if n_Q <= n_T) + + for n_Q in n_Q_values: + Q = np.random.rand(n_Q) + T = np.random.rand(n_T) + ref = naive.rolling_window_dot_product(Q, T) + for func_name in get_sdp_function_names(): + func = getattr(sdp, func_name) + try: + comp = func(Q, T) + npt.assert_allclose(comp, ref) + except Exception as e: # pragma: no cover + msg = f"Error in {func_name}, with n_Q={len(Q)} and n_T={len(T)}" + warnings.warn(msg) + raise e + + +def test_sdp_power2(): + # test for case 5. len(T) is power of 2 + pmin = 3 + pmax = 13 + for func_name in get_sdp_function_names(): func = getattr(sdp, func_name) try: - comp = func(Q, T) - ref = naive.rolling_window_dot_product(Q, T) - npt.assert_allclose(comp, ref) + for q in range(pmin, pmax + 1): + n_Q = 2**q + for p in range(q, pmax + 1): + n_T = 2**p + Q = np.random.rand(n_Q) + T = np.random.rand(n_T) + + ref = naive.rolling_window_dot_product(Q, T) + comp = func(Q, T) + npt.assert_allclose(comp, ref) + except Exception as e: # pragma: no cover - msg = f"Error in {func_name}, with n_Q={len(Q)} and n_T={len(T)}" + msg = f"Error in {func_name}, with q={q} and p={p}" warnings.warn(msg) raise e + + return From 275a08786c1a5e0781d2f3a67a52dcfc8c0c276d Mon Sep 17 00:00:00 2001 From: NimaSarajpoor Date: Wed, 11 Feb 2026 02:03:26 -0500 Subject: [PATCH 10/11] add minor comment --- tests/test_sdp.py | 1 + 1 file changed, 1 insertion(+) diff --git a/tests/test_sdp.py b/tests/test_sdp.py index 3f9b18949..53aebac1a 100644 --- a/tests/test_sdp.py +++ b/tests/test_sdp.py @@ -37,6 +37,7 @@ # 1000-2000. This allows us to test sliding_dot_product for # wider range of query lengths. +# test cases 1-4 test_inputs = [ # Input format: # ( From 55589d3fe55862d4de901fa8fb5a803b1ac5fa71 Mon Sep 17 00:00:00 2001 From: Nima Sarajpoor Date: Wed, 11 Feb 2026 08:43:02 -0500 Subject: [PATCH 11/11] Minor change --- tests/test_sdp.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/tests/test_sdp.py b/tests/test_sdp.py index 53aebac1a..b8d98eeb9 100644 --- a/tests/test_sdp.py +++ b/tests/test_sdp.py @@ -122,7 +122,7 @@ def test_sdp(n_T, remainder, comparator): comp = func(Q, T) npt.assert_allclose(comp, ref) except Exception as e: # pragma: no cover - msg = f"Error in {func_name}, with n_Q={len(Q)} and n_T={len(T)}" + msg = f"Error in {func_name}, with n_Q={n_Q} and n_T={n_T}" warnings.warn(msg) raise e