@@ -49,54 +49,77 @@ namespace deglib::distances::fp32_ip {
4949
5050#if defined(DEGLIB_X86)
5151 // -------------------------------------------------------------------
52- // InnerProductFloat16Ext — processes 16 floats (64 bytes) per iteration.
52+ // InnerProductFloat SIMD implementations — process vectors with
53+ // aligned SIMD portions plus scalar residuals for any unaligned tail.
5354 // Separate classes per SIMD width so that compare() has zero
5455 // runtime dispatch overhead — select_dist() chooses the class.
56+ // The HasResidual template parameter controls whether the scalar
57+ // residual tail loop is compiled in. When HasResidual == false,
58+ // the residual loop is eliminated at compile time, producing a
59+ // faster path for dimensions that are known to be SIMD-aligned.
5560 // -------------------------------------------------------------------
5661
57- class InnerProductFloat16Ext_AVX512 {
62+ template <ResidualMode Mode = ResidualMode::Full>
63+ class InnerProductFloat_AVX512 {
64+ static constexpr bool HasDualSimd = has_flag(Mode, ResidualMode::DualSimd);
65+ static constexpr bool HasSimd = has_flag(Mode, ResidualMode::Simd);
66+ static constexpr bool HasTail = has_flag(Mode, ResidualMode::Tail);
67+
5868 public:
5969 DEGLIB_TARGET_AVX512 inline static float compare (const void *pVect1v, const void *pVect2v, const void *qty_ptr) {
60- return 1 .f - dot (pVect1v, pVect2v, qty_ptr);
61- }
62-
63- DEGLIB_TARGET_AVX512 inline static float dot (const void *pVect1v, const void *pVect2v, const void *qty_ptr) {
64- return ip_16ext (pVect1v, pVect2v, qty_ptr);
65- }
66-
67- DEGLIB_TARGET_AVX512 inline static float ip_16ext (const void *pVect1v, const void *pVect2v, const void *qty_ptr) {
6870 float *a = (float *) pVect1v;
6971 float *b = (float *) pVect2v;
7072 size_t size = *((size_t *) qty_ptr);
7173
7274 const float *last = a + size;
7375
74- __m512 sum512 = _mm512_setzero_ps ();
75- while (a < last) {
76- sum512 = _mm512_fmadd_ps (_mm512_loadu_ps (a), _mm512_loadu_ps (b), sum512);
77- a += 16 ;
78- b += 16 ;
76+ __m512 sum512_1 = _mm512_setzero_ps ();
77+ __m512 sum512_2 = _mm512_setzero_ps ();
78+ if constexpr (HasDualSimd) {
79+ while (a + 31 < last) {
80+ sum512_1 = _mm512_fmadd_ps (_mm512_loadu_ps (a), _mm512_loadu_ps (b), sum512_1);
81+ a += 16 ;
82+ b += 16 ;
83+ sum512_2 = _mm512_fmadd_ps (_mm512_loadu_ps (a), _mm512_loadu_ps (b), sum512_2);
84+ a += 16 ;
85+ b += 16 ;
86+ }
87+ }
88+ if constexpr (HasSimd) {
89+ while (a + 15 < last) {
90+ sum512_1 = _mm512_fmadd_ps (_mm512_loadu_ps (a), _mm512_loadu_ps (b), sum512_1);
91+ a += 16 ;
92+ b += 16 ;
93+ }
7994 }
8095
96+ // Horizontal reduce of SIMD accumulators
97+ __m512 sum512 = _mm512_add_ps (sum512_1, sum512_2);
8198 __m256 sum256 = _mm256_add_ps (_mm512_extractf32x8_ps (sum512, 0 ), _mm512_extractf32x8_ps (sum512, 1 ));
8299 __m128 sum128 = _mm_add_ps (_mm256_extractf128_ps (sum256, 0 ), _mm256_extractf128_ps (sum256, 1 ));
83100 alignas (32 ) float f[4 ];
84101 _mm_store_ps (f, sum128);
85- return f[0 ] + f[1 ] + f[2 ] + f[3 ];
102+ float result = f[0 ] + f[1 ] + f[2 ] + f[3 ];
103+
104+ // Scalar residual for the unaligned tail — eliminated at compile-time if HasTail == false
105+ if constexpr (HasTail) {
106+ while (a < last) {
107+ result = std::fma (*a++, *b++, result);
108+ }
109+ }
110+
111+ return 1 .f - result;
86112 }
87113 };
88114
89- class InnerProductFloat16Ext_AVX2 {
115+ template <ResidualMode Mode = ResidualMode::Full>
116+ class InnerProductFloat_AVX2 {
117+ static constexpr bool HasDualSimd = has_flag(Mode, ResidualMode::DualSimd);
118+ static constexpr bool HasSimd = has_flag(Mode, ResidualMode::Simd);
119+ static constexpr bool HasTail = has_flag(Mode, ResidualMode::Tail);
120+
90121 public:
91122 DEGLIB_TARGET_AVX2 inline static float compare (const void *pVect1v, const void *pVect2v, const void *qty_ptr) {
92- return 1 .f - dot (pVect1v, pVect2v, qty_ptr);
93- }
94-
95- DEGLIB_TARGET_AVX2 inline static float dot (const void *pVect1v, const void *pVect2v, const void *qty_ptr) {
96- return ip_16ext (pVect1v, pVect2v, qty_ptr);
97- }
98-
99- DEGLIB_TARGET_AVX2 inline static float ip_16ext (const void *pVect1v, const void *pVect2v, const void *qty_ptr) {
100123 float *a = (float *) pVect1v;
101124 float *b = (float *) pVect2v;
102125 size_t size = *((size_t *) qty_ptr);
@@ -105,131 +128,88 @@ namespace deglib::distances::fp32_ip {
105128
106129 __m256 sum256_1 = _mm256_setzero_ps ();
107130 __m256 sum256_2 = _mm256_setzero_ps ();
108- while (a < last) {
109- sum256_1 = _mm256_fmadd_ps (_mm256_loadu_ps (a), _mm256_loadu_ps (b), sum256_1);
110- a += 8 ;
111- b += 8 ;
112- sum256_2 = _mm256_fmadd_ps (_mm256_loadu_ps (a), _mm256_loadu_ps (b), sum256_2);
113- a += 8 ;
114- b += 8 ;
131+ if constexpr (HasDualSimd) {
132+ while (a + 15 < last) {
133+ sum256_1 = _mm256_fmadd_ps (_mm256_loadu_ps (a), _mm256_loadu_ps (b), sum256_1);
134+ a += 8 ;
135+ b += 8 ;
136+ sum256_2 = _mm256_fmadd_ps (_mm256_loadu_ps (a), _mm256_loadu_ps (b), sum256_2);
137+ a += 8 ;
138+ b += 8 ;
139+ }
140+ }
141+ if constexpr (HasSimd) {
142+ while (a + 7 < last) {
143+ sum256_1 = _mm256_fmadd_ps (_mm256_loadu_ps (a), _mm256_loadu_ps (b), sum256_1);
144+ a += 8 ;
145+ b += 8 ;
146+ }
115147 }
148+
149+ // Horizontal reduce of SIMD accumulators
116150 __m256 sum256 = _mm256_add_ps (sum256_1, sum256_2);
117151 __m128 sum128 = _mm_add_ps (_mm256_extractf128_ps (sum256, 0 ), _mm256_extractf128_ps (sum256, 1 ));
118152 alignas (32 ) float f[4 ];
119153 _mm_store_ps (f, sum128);
120- return f[0 ] + f[1 ] + f[2 ] + f[3 ];
121- }
122- };
123-
124-
125- // -------------------------------------------------------------------
126- // InnerProductFloat8Ext — processes 8 floats (32 bytes) per iteration.
127- // -------------------------------------------------------------------
154+ float result = f[0 ] + f[1 ] + f[2 ] + f[3 ];
128155
129- class InnerProductFloat8Ext_AVX2 {
130- public:
131- DEGLIB_TARGET_AVX2 inline static float compare (const void *pVect1v, const void *pVect2v, const void *qty_ptr) {
132- return 1 .f - dot (pVect1v, pVect2v, qty_ptr);
133- }
134-
135- DEGLIB_TARGET_AVX2 inline static float dot (const void *pVect1v, const void *pVect2v, const void *qty_ptr) {
136- return ip_8ext (pVect1v, pVect2v, qty_ptr);
137- }
138-
139- DEGLIB_TARGET_AVX2 inline static float ip_8ext (const void *pVect1v, const void *pVect2v, const void *qty_ptr) {
140- float *a = (float *) pVect1v;
141- float *b = (float *) pVect2v;
142- size_t size = *((size_t *) qty_ptr);
143-
144- const float *last = a + size;
145-
146- __m256 sum256 = _mm256_setzero_ps ();
147- while (a < last) {
148- sum256 = _mm256_fmadd_ps (_mm256_loadu_ps (a), _mm256_loadu_ps (b), sum256);
149- a += 8 ;
150- b += 8 ;
156+ // Scalar residual for the unaligned tail — eliminated at compile-time if HasTail == false
157+ if constexpr (HasTail) {
158+ while (a < last) {
159+ result = std::fma (*a++, *b++, result);
160+ }
151161 }
152- __m128 sum128 = _mm_add_ps (_mm256_extractf128_ps (sum256, 0 ), _mm256_extractf128_ps (sum256, 1 ));
153- alignas (32 ) float f[4 ];
154- _mm_store_ps (f, sum128);
155- return f[0 ] + f[1 ] + f[2 ] + f[3 ];
156- }
157- };
158-
159-
160- // -------------------------------------------------------------------
161- // InnerProductFloat4Ext — processes 4 floats (16 bytes) per iteration.
162- // -------------------------------------------------------------------
163-
164162
165- // -------------------------------------------------------------------
166- // Residual classes — process the aligned portion with the SIMD
167- // variant and the remainder with the scalar fallback.
168- // -------------------------------------------------------------------
169-
170- class InnerProductFloat16ExtResiduals_AVX512 {
171- public:
172- DEGLIB_TARGET_AVX512 inline static float compare (const void *pVect1v, const void *pVect2v, const void *qty_ptr) {
173- size_t qty = *((size_t *) qty_ptr);
174-
175- size_t qty16 = qty >> 4 << 4 ;
176- float res = InnerProductFloat16Ext_AVX512::dot (pVect1v, pVect2v, &qty16);
177- float *pVect1 = (float *) pVect1v + qty16;
178- float *pVect2 = (float *) pVect2v + qty16;
179-
180- size_t qty_left = qty - qty16;
181- float res_tail = InnerProductFloat::dot (pVect1, pVect2, &qty_left);
182- return 1 .f - (res + res_tail);
183- }
184- };
185-
186- class InnerProductFloat16ExtResiduals_AVX2 {
187- public:
188- DEGLIB_TARGET_AVX2 inline static float compare (const void *pVect1v, const void *pVect2v, const void *qty_ptr) {
189- size_t qty = *((size_t *) qty_ptr);
190-
191- size_t qty16 = qty >> 4 << 4 ;
192- float res = InnerProductFloat16Ext_AVX2::dot (pVect1v, pVect2v, &qty16);
193- float *pVect1 = (float *) pVect1v + qty16;
194- float *pVect2 = (float *) pVect2v + qty16;
195-
196- size_t qty_left = qty - qty16;
197- float res_tail = InnerProductFloat::dot (pVect1, pVect2, &qty_left);
198- return 1 .f - (res + res_tail);
163+ return 1 .f - result;
199164 }
200165 };
201-
202-
203166#endif
204167
205168 using DistanceVariant = std::variant<
206169 InnerProductFloat
207170#if defined(DEGLIB_X86)
208- , InnerProductFloat16Ext_AVX512,
209- InnerProductFloat16Ext_AVX2,
210- InnerProductFloat8Ext_AVX2,
211- InnerProductFloat16ExtResiduals_AVX512,
212- InnerProductFloat16ExtResiduals_AVX2
171+ , InnerProductFloat_AVX512<ResidualMode::Full>
172+ , InnerProductFloat_AVX512<ResidualMode::DualPlusSimd>
173+ , InnerProductFloat_AVX512<ResidualMode::DualTail>
174+ , InnerProductFloat_AVX512<ResidualMode::DualOnly>
175+ , InnerProductFloat_AVX512<ResidualMode::SimdTail>
176+ , InnerProductFloat_AVX512<ResidualMode::SimdOnly>
177+ , InnerProductFloat_AVX512<ResidualMode::TailOnly>
178+ , InnerProductFloat_AVX2<ResidualMode::Full>
179+ , InnerProductFloat_AVX2<ResidualMode::DualPlusSimd>
180+ , InnerProductFloat_AVX2<ResidualMode::DualTail>
181+ , InnerProductFloat_AVX2<ResidualMode::DualOnly>
182+ , InnerProductFloat_AVX2<ResidualMode::SimdTail>
183+ , InnerProductFloat_AVX2<ResidualMode::SimdOnly>
184+ , InnerProductFloat_AVX2<ResidualMode::TailOnly>
213185#endif
214186 >;
215187
216188 inline DistanceVariant select_dist (const size_t dim) {
217189#if defined(DEGLIB_X86)
218190 if (deglib::cpu::has_avx512 ()) {
219- if (dim % 16 == 0 )
220- return InnerProductFloat16Ext_AVX512{};
221- else if (dim % 8 == 0 )
222- return InnerProductFloat8Ext_AVX2{};
223- else if (dim > 16 )
224- return InnerProductFloat16ExtResiduals_AVX512{};
191+ if (dim < 16 ) {
192+ return InnerProductFloat_AVX512<ResidualMode::TailOnly>{};
193+ } else if (dim < 32 ) {
194+ if (dim == 16 ) return InnerProductFloat_AVX512<ResidualMode::SimdOnly>{};
195+ else return InnerProductFloat_AVX512<ResidualMode::SimdTail>{};
196+ } else {
197+ if (dim % 32 == 0 ) return InnerProductFloat_AVX512<ResidualMode::DualOnly>{};
198+ else if (dim % 16 == 0 ) return InnerProductFloat_AVX512<ResidualMode::DualPlusSimd>{};
199+ else return InnerProductFloat_AVX512<ResidualMode::Full>{};
200+ }
225201 }
226202 else if (deglib::cpu::has_avx2 ()) {
227- if (dim % 16 == 0 )
228- return InnerProductFloat16Ext_AVX2{};
229- else if (dim % 8 == 0 )
230- return InnerProductFloat8Ext_AVX2{};
231- else if (dim > 16 )
232- return InnerProductFloat16ExtResiduals_AVX2{};
203+ if (dim < 8 ) {
204+ return InnerProductFloat_AVX2<ResidualMode::TailOnly>{};
205+ } else if (dim < 16 ) {
206+ if (dim == 8 ) return InnerProductFloat_AVX2<ResidualMode::SimdOnly>{};
207+ else return InnerProductFloat_AVX2<ResidualMode::SimdTail>{};
208+ } else {
209+ if (dim % 16 == 0 ) return InnerProductFloat_AVX2<ResidualMode::DualOnly>{};
210+ else if (dim % 8 == 0 ) return InnerProductFloat_AVX2<ResidualMode::DualPlusSimd>{};
211+ else return InnerProductFloat_AVX2<ResidualMode::Full>{};
212+ }
233213 }
234214#endif
235215 return InnerProductFloat{};
0 commit comments