ml: SIMD for SVM kernel reductions

SVM::predict is dominated by the per-feature kernel reduction over the
support vectors. Vectorize the four reduction kernels in svm.cpp with
universal intrinsics (two independent accumulators + scalar tail):
calc_non_rbf_base (dot product), calc_rbf (squared distance),
calc_intersec (min-sum) and calc_chi2. Same approach as the KNN
findNearest reduction in #29380. Adds modules/ml/perf/perf_svm.cpp
covering the RBF/POLY/INTER/CHI2 kernels.

modules/ml/perf/perf_svm.cpp

@@ -0,0 +1,48 @@
+// This file is part of OpenCV project.
+// It is subject to the license terms in the LICENSE file found in the top-level directory
+// of this distribution and at http://opencv.org/license.html.
+#include "perf_precomp.hpp"
+
+namespace opencv_test { namespace {
+
+// SVM::predict is dominated by the per-feature kernel reduction over the support
+// vectors: calc_non_rbf_base (LINEAR/POLY/SIGMOID dot product), calc_rbf
+// (squared distance), calc_intersec (min-sum) and calc_chi2. Train data is kept
+// non-negative so the INTER and CHI2 kernels stay in their valid domain.
+typedef TestBaseWithParam< tuple<int, int, int> > SVMPredict;  // (samples, dims, kernelType)
+
+PERF_TEST_P(SVMPredict, kernels, testing::Combine(
+            testing::Values(1500),
+            testing::Values(512),
+            testing::Values((int)SVM::RBF, (int)SVM::POLY, (int)SVM::INTER, (int)SVM::CHI2)))
+{
+    const int nsamples = get<0>(GetParam());
+    const int dims     = get<1>(GetParam());
+    const int kernel   = get<2>(GetParam());
+    const int nquery   = 1000;
+
+    Mat train(nsamples, dims, CV_32F), query(nquery, dims, CV_32F);
+    Mat responses(nsamples, 1, CV_32S);
+    RNG& rng = theRNG();
+    rng.fill(train, RNG::UNIFORM, 0.f, 1.f);
+    rng.fill(query, RNG::UNIFORM, 0.f, 1.f);
+    for (int i = 0; i < nsamples; i++)
+        responses.at<int>(i) = i & 1;
+
+    Ptr<SVM> svm = SVM::create();
+    svm->setType(SVM::C_SVC);
+    svm->setKernel(kernel);
+    svm->setC(1);
+    svm->setGamma(0.1);
+    svm->setDegree(3);   // POLY
+    svm->setCoef0(1);    // POLY
+    svm->setTermCriteria(TermCriteria(TermCriteria::MAX_ITER, 1000, 1e-3));
+    svm->train(train, ROW_SAMPLE, responses);
+
+    Mat results;
+    TEST_CYCLE() svm->predict(query, results);
+
+    SANITY_CHECK_NOTHING();
+}
+
+}} // namespace

modules/ml/src/svm.cpp

@@ -41,6 +41,7 @@
 //M*/
 
 #include "precomp.hpp"
+#include "opencv2/core/hal/intrin.hpp"
 #include <opencv2/core/utils/logger.hpp>
 
 #include <stdarg.h>
@@ -173,9 +174,19 @@ public:
         {
             const float* sample = &vecs[j*var_count];
             double s = 0;
-            for( k = 0; k <= var_count - 4; k += 4 )
+            k = 0;
-                s += sample[k]*another[k] + sample[k+1]*another[k+1] +
+#if (CV_SIMD || CV_SIMD_SCALABLE)
-                sample[k+2]*another[k+2] + sample[k+3]*another[k+3];
+            {
+                const int vl = VTraits<v_float32>::vlanes();
+                v_float32 a0 = vx_setzero_f32(), a1 = vx_setzero_f32();
+                for( ; k <= var_count - 2*vl; k += 2*vl )
+                {
+                    a0 = v_fma(vx_load(sample + k), vx_load(another + k), a0);
+                    a1 = v_fma(vx_load(sample + k + vl), vx_load(another + k + vl), a1);
+                }
+                s = (double)v_reduce_sum(v_add(a0, a1));
+            }
+#endif
             for( ; k < var_count; k++ )
                 s += sample[k]*another[k];
             results[j] = (Qfloat)(s*alpha + beta);
@@ -228,20 +239,21 @@ public:
         {
             const float* sample = &vecs[j*var_count];
             double s = 0;
-
+            k = 0;
-            for( k = 0; k <= var_count - 4; k += 4 )
+#if (CV_SIMD || CV_SIMD_SCALABLE)
             {
-                double t0 = sample[k] - another[k];
+                const int vl = VTraits<v_float32>::vlanes();
-                double t1 = sample[k+1] - another[k+1];
+                v_float32 a0 = vx_setzero_f32(), a1 = vx_setzero_f32();
-
+                for( ; k <= var_count - 2*vl; k += 2*vl )
-                s += t0*t0 + t1*t1;
+                {
-
+                    v_float32 d0 = v_sub(vx_load(sample + k), vx_load(another + k));
-                t0 = sample[k+2] - another[k+2];
+                    a0 = v_fma(d0, d0, a0);
-                t1 = sample[k+3] - another[k+3];
+                    v_float32 d1 = v_sub(vx_load(sample + k + vl), vx_load(another + k + vl));
-
+                    a1 = v_fma(d1, d1, a1);
-                s += t0*t0 + t1*t1;
+                }
+                s = (double)v_reduce_sum(v_add(a0, a1));
             }
-
+#endif
             for( ; k < var_count; k++ )
             {
                 double t0 = sample[k] - another[k];
@@ -266,9 +278,19 @@ public:
         {
             const float* sample = &vecs[j*var_count];
             double s = 0;
-            for( k = 0; k <= var_count - 4; k += 4 )
+            k = 0;
-                s += std::min(sample[k],another[k]) + std::min(sample[k+1],another[k+1]) +
+#if (CV_SIMD || CV_SIMD_SCALABLE)
-                std::min(sample[k+2],another[k+2]) + std::min(sample[k+3],another[k+3]);
+            {
+                const int vl = VTraits<v_float32>::vlanes();
+                v_float32 a0 = vx_setzero_f32(), a1 = vx_setzero_f32();
+                for( ; k <= var_count - 2*vl; k += 2*vl )
+                {
+                    a0 = v_add(a0, v_min(vx_load(sample + k), vx_load(another + k)));
+                    a1 = v_add(a1, v_min(vx_load(sample + k + vl), vx_load(another + k + vl)));
+                }
+                s = (double)v_reduce_sum(v_add(a0, a1));
+            }
+#endif
             for( ; k < var_count; k++ )
                 s += std::min(sample[k],another[k]);
             results[j] = (Qfloat)(s);
@@ -286,7 +308,23 @@ public:
         {
             const float* sample = &vecs[j*var_count];
             double chi2 = 0;
-            for(k = 0 ; k < var_count; k++ )
+            k = 0;
+#if (CV_SIMD || CV_SIMD_SCALABLE)
+            {
+                const int vl = VTraits<v_float32>::vlanes();
+                v_float32 acc = vx_setzero_f32(), z = vx_setzero_f32();
+                for( ; k <= var_count - vl; k += vl )
+                {
+                    v_float32 a = vx_load(sample + k), b = vx_load(another + k);
+                    v_float32 d = v_sub(a, b), dv = v_add(a, b);
+                    v_float32 term = v_div(v_mul(d, d), dv);
+                    term = v_select(v_eq(dv, z), z, term);  // divisor==0 -> add 0 (skip)
+                    acc = v_add(acc, term);
+                }
+                chi2 = (double)v_reduce_sum(acc);
+            }
+#endif
+            for( ; k < var_count; k++ )
             {
                 double d = sample[k]-another[k];
                 double devisor = sample[k]+another[k];