char-fork/llama.cpp

git clone https://git.t4t.associates/char-fork/llama.cpp

cqderekggml-hexagon: add F16 support for unary ops (#28228)f027c4f1b

char/aliceai-t5-moe
19.4 KiB454 linesraw
1#ifndef HVX_NORM_H
2#define HVX_NORM_H
3
4#include <stdint.h>
5#include "hvx-base.h"
6#include "hvx-reduce.h"
7#include "hvx-inverse.h"
8#include "hvx-sqrt.h"
9#include "hvx-repl.h"
10
11static inline void hvx_fast_rms_norm_f32(const uint8_t * restrict src,
12                                         uint8_t * restrict dst,
13                                         const int num_elems,
14                                         float     epsilon) {
15
16    const HVX_Vector * restrict v_src = (HVX_Vector *) src;
17    HVX_Vector * restrict v_dst       = (HVX_Vector *) dst;
18
19    const int nvec = num_elems / VLEN_FP32;    // number of full vectors
20    const int nloe = num_elems % VLEN_FP32;    // leftover elements
21
22    // Compute sum of squares for full vectors
23    HVX_Vector sum_v = Q6_V_vsplat_R(0x00000000);
24    HVX_Vector epsilon_v = hvx_vec_splat_f32(epsilon);
25
26    #pragma unroll(4)
27    for (int i = 0; i < nvec; i++) {
28        HVX_Vector v1 = v_src[i];
29        HVX_Vector v2 = Q6_Vqf32_vmpy_VsfVsf(v1, v1);
30        sum_v = Q6_Vqf32_vadd_Vqf32Vqf32(sum_v, v2);
31    }
32
33    // Handle tail elements using vectorized ops with masking
34    if (nloe > 0) {
35        HVX_VectorPred bmask = Q6_Q_vsetq_R(nloe * 4);
36        HVX_Vector v1 = Q6_V_vand_QV(bmask, v_src[nvec]);
37        HVX_Vector v2 = Q6_Vqf32_vmpy_VsfVsf(v1, v1);
38        sum_v = Q6_Vqf32_vadd_Vqf32Vqf32(sum_v, v2);
39    }
40
41    // Reduce HVX sum
42    sum_v = hvx_vec_reduce_sum_f32(Q6_Vsf_equals_Vqf32(sum_v));
43
44    HVX_Vector t_v            = hvx_vec_splat_f32((float) num_elems);
45    HVX_Vector denom_v        = hvx_vec_inverse_f32(t_v);
46    HVX_Vector mean_v         = Q6_Vqf32_vmpy_VsfVsf(sum_v, denom_v);
47    HVX_Vector mean_epsilon_v = Q6_Vqf32_vadd_Vqf32Vsf(mean_v, epsilon_v);
48
49    // Scale full vectors
50    HVX_Vector scale_v = hvx_vec_rsqrt_f32(Q6_Vsf_equals_Vqf32(mean_epsilon_v));
51
52    #pragma unroll(4)
53    for (int i = 0; i < nvec; i++) {
54        HVX_Vector v1 = v_src[i];
55        HVX_Vector v2 = Q6_Vqf32_vmpy_VsfVsf(v1, scale_v);
56        v_dst[i] = Q6_Vsf_equals_Vqf32(v2);
57    }
58
59    // Handle tail elements using vectorized ops with masking
60    if (nloe > 0) {
61        HVX_VectorPred bmask = Q6_Q_vsetq_R(nloe * 4);
62        HVX_Vector v1 = Q6_V_vand_QV(bmask, v_src[nvec]);
63        HVX_Vector v2 = Q6_Vqf32_vmpy_VsfVsf(v1, scale_v);
64        HVX_Vector result = Q6_Vsf_equals_Vqf32(v2);
65
66        // Store with masking to avoid overwriting memory beyond the tensor
67        hvx_vec_store_a(&v_dst[nvec], nloe * 4, result);
68    }
69}
70
71static inline void hvx_fast_rms_norm_mul_f32(const uint8_t * restrict src,
72                                             const uint8_t * restrict weight,
73                                             uint8_t * restrict dst,
74                                             const int num_elems,
75                                             float     epsilon) {
76    const HVX_Vector * restrict v_src    = (const HVX_Vector *) src;
77    const HVX_Vector * restrict v_weight = (const HVX_Vector *) weight;
78    HVX_Vector * restrict v_dst          = (HVX_Vector *) dst;
79
80    const int nvec = num_elems / VLEN_FP32;    // number of full vectors
81    const int nloe = num_elems % VLEN_FP32;    // leftover elements
82
83    // Compute sum of squares for full vectors
84    HVX_Vector sum_v = Q6_V_vsplat_R(0x00000000);
85    HVX_Vector epsilon_v = hvx_vec_splat_f32(epsilon);
86
87    #pragma unroll(4)
88    for (int i = 0; i < nvec; i++) {
89        HVX_Vector v1 = v_src[i];
90        HVX_Vector v2 = Q6_Vqf32_vmpy_VsfVsf(v1, v1);
91        sum_v = Q6_Vqf32_vadd_Vqf32Vqf32(sum_v, v2);
92    }
93
94    // Handle tail elements using vectorized ops with masking
95    if (nloe > 0) {
96        HVX_VectorPred bmask = Q6_Q_vsetq_R(nloe * 4);
97        HVX_Vector v1 = Q6_V_vand_QV(bmask, v_src[nvec]);
98        HVX_Vector v2 = Q6_Vqf32_vmpy_VsfVsf(v1, v1);
99        sum_v = Q6_Vqf32_vadd_Vqf32Vqf32(sum_v, v2);
100    }
101
102    // Reduce HVX sum
103    sum_v = hvx_vec_reduce_sum_f32(Q6_Vsf_equals_Vqf32(sum_v));
104
105    HVX_Vector t_v            = hvx_vec_splat_f32((float) num_elems);
106    HVX_Vector denom_v        = hvx_vec_inverse_f32(t_v);
107    HVX_Vector mean_v         = Q6_Vqf32_vmpy_VsfVsf(sum_v, denom_v);
108    HVX_Vector mean_epsilon_v = Q6_Vqf32_vadd_Vqf32Vsf(mean_v, epsilon_v);
109
110    // Scale and multiply
111    HVX_Vector scale_v = hvx_vec_rsqrt_f32(Q6_Vsf_equals_Vqf32(mean_epsilon_v));
112
113    #pragma unroll(4)
114    for (int i = 0; i < nvec; i++) {
115        HVX_Vector v1 = v_src[i];
116        HVX_Vector v2 = Q6_Vqf32_vmpy_VsfVsf(v1, scale_v);
117        HVX_Vector v3 = Q6_Vsf_equals_Vqf32(v2);
118        HVX_Vector result = Q6_Vqf32_vmpy_VsfVsf(v3, v_weight[i]);
119        v_dst[i] = Q6_Vsf_equals_Vqf32(result);
120    }
121
122    // Handle tail elements using vectorized ops with masking
123    if (nloe > 0) {
124        HVX_VectorPred bmask = Q6_Q_vsetq_R(nloe * 4);
125        HVX_Vector v1 = Q6_V_vand_QV(bmask, v_src[nvec]);
126        HVX_Vector v2 = Q6_Vqf32_vmpy_VsfVsf(v1, scale_v);
127        HVX_Vector v3 = Q6_Vsf_equals_Vqf32(v2);
128        HVX_Vector result = Q6_Vqf32_vmpy_VsfVsf(v3, v_weight[nvec]);
129        HVX_Vector res_v = Q6_Vsf_equals_Vqf32(result);
130
131        // Store with masking to avoid overwriting memory beyond the tensor
132        hvx_vec_store_a(&v_dst[nvec], nloe * 4, res_v);
133    }
134}
135
136static inline void hvx_fast_norm_f32(const uint8_t * restrict src,
137                                     uint8_t * restrict dst,
138                                     const int num_elems,
139                                     float     epsilon) {
140
141    const HVX_Vector * restrict v_src = (HVX_Vector *) src;
142    HVX_Vector * restrict v_dst       = (HVX_Vector *) dst;
143
144    const int nvec = num_elems / VLEN_FP32;    // number of full vectors
145    const int nloe = num_elems % VLEN_FP32;    // leftover elements
146
147    // Compute sum of squares and sum of values for full vectors
148    HVX_Vector sum_sq_v = Q6_V_vsplat_R(0x00000000);
149    HVX_Vector sum_x_v  = Q6_V_vsplat_R(0x00000000);
150    HVX_Vector epsilon_v = hvx_vec_splat_f32(epsilon);
151
152    #pragma unroll(4)
153    for (int i = 0; i < nvec; i++) {
154        HVX_Vector v1 = v_src[i];
155        HVX_Vector v2 = Q6_Vqf32_vmpy_VsfVsf(v1, v1);
156        sum_sq_v = Q6_Vqf32_vadd_Vqf32Vqf32(sum_sq_v, v2);
157        sum_x_v  = Q6_Vqf32_vadd_Vqf32Vqf32(sum_x_v,  Q6_Vqf32_vadd_VsfVsf(v1, Q6_V_vzero()));
158    }
159
160    // Handle tail elements using vectorized ops with masking
161    if (nloe > 0) {
162        HVX_VectorPred bmask = Q6_Q_vsetq_R(nloe * 4);
163        HVX_Vector v1 = Q6_V_vand_QV(bmask, v_src[nvec]);
164        HVX_Vector v2 = Q6_Vqf32_vmpy_VsfVsf(v1, v1);
165        sum_sq_v = Q6_Vqf32_vadd_Vqf32Vqf32(sum_sq_v, v2);
166        sum_x_v  = Q6_Vqf32_vadd_Vqf32Vqf32(sum_x_v,  Q6_Vqf32_vadd_VsfVsf(v1, Q6_V_vzero()));
167    }
168
169    // Reduce HVX sums
170    sum_sq_v = hvx_vec_reduce_sum_f32(Q6_Vsf_equals_Vqf32(sum_sq_v));
171    sum_x_v  = hvx_vec_reduce_sum_f32(Q6_Vsf_equals_Vqf32(sum_x_v));
172
173    HVX_Vector t_v            = hvx_vec_splat_f32((float) num_elems);
174    HVX_Vector denom_v        = hvx_vec_inverse_f32(t_v);
175    HVX_Vector mean_sq_v      = Q6_Vqf32_vmpy_VsfVsf(sum_sq_v, denom_v);
176    HVX_Vector mean_x_v       = Q6_Vqf32_vmpy_VsfVsf(sum_x_v,  denom_v);
177    HVX_Vector mean_x_sq_v    = Q6_Vqf32_vmpy_VsfVsf(Q6_Vsf_equals_Vqf32(mean_x_v), Q6_Vsf_equals_Vqf32(mean_x_v));
178    HVX_Vector var_v          = Q6_Vqf32_vsub_Vqf32Vqf32(mean_sq_v, mean_x_sq_v);
179    HVX_Vector var_epsilon_v  = Q6_Vqf32_vadd_Vqf32Vsf(var_v, epsilon_v);
180
181    // scale = rsqrt(variance + epsilon),  mean_x broadcast for subtraction
182    HVX_Vector scale_v  = hvx_vec_rsqrt_f32(Q6_Vsf_equals_Vqf32(var_epsilon_v));
183    HVX_Vector mean_x_b = hvx_vec_repl_f32(Q6_Vsf_equals_Vqf32(mean_x_v));
184
185    #pragma unroll(4)
186    for (int i = 0; i < nvec; i++) {
187        HVX_Vector v1 = v_src[i];
188        HVX_Vector v2 = Q6_Vqf32_vsub_VsfVsf(v1, mean_x_b);
189        HVX_Vector v3 = Q6_Vqf32_vmpy_VsfVsf(Q6_Vsf_equals_Vqf32(v2), scale_v);
190        v_dst[i] = Q6_Vsf_equals_Vqf32(v3);
191    }
192
193    // Handle tail elements using vectorized ops with masking
194    if (nloe > 0) {
195        HVX_VectorPred bmask = Q6_Q_vsetq_R(nloe * 4);
196        HVX_Vector v1 = Q6_V_vand_QV(bmask, v_src[nvec]);
197        HVX_Vector v2 = Q6_Vqf32_vsub_VsfVsf(v1, mean_x_b);
198        HVX_Vector v3 = Q6_Vqf32_vmpy_VsfVsf(Q6_Vsf_equals_Vqf32(v2), scale_v);
199        HVX_Vector result = Q6_Vsf_equals_Vqf32(v3);
200
201        // Store with masking to avoid overwriting memory beyond the tensor
202        hvx_vec_store_a(&v_dst[nvec], nloe * 4, result);
203    }
204}
205
206static inline void hvx_fast_l2_norm_f32(const uint8_t * restrict src,
207                                        uint8_t * restrict dst,
208                                        const int num_elems,
209                                        float     epsilon) {
210
211    const HVX_Vector * restrict v_src = (HVX_Vector *) src;
212    HVX_Vector * restrict v_dst       = (HVX_Vector *) dst;
213
214    HVX_Vector sum_v = hvx_vec_splat_f32(0.0f);
215
216    const int nvec = num_elems / VLEN_FP32;
217    const int nloe = num_elems % VLEN_FP32;
218
219    #pragma unroll(4)
220    for (int i = 0; i < nvec; i++) {
221        HVX_Vector v1 = v_src[i];
222        HVX_Vector sq = Q6_Vqf32_vmpy_VsfVsf(v1, v1);
223        sum_v         = Q6_Vqf32_vadd_Vqf32Vqf32(sum_v, sq);
224    }
225
226    // Include tail elements in the sum-of-squares using a predicate mask
227    if (nloe > 0) {
228        HVX_VectorPred bmask = Q6_Q_vsetq_R(nloe * 4);
229        HVX_Vector v1 = Q6_V_vand_QV(bmask, v_src[nvec]);
230        HVX_Vector sq = Q6_Vqf32_vmpy_VsfVsf(v1, v1);
231        sum_v         = Q6_Vqf32_vadd_Vqf32Vqf32(sum_v, sq);
232    }
233
234    // Compute scale = 1/fmax(sqrt(sum), epsilon) entirely in HVX registers.
235    // hvx_vec_rsqrt_f32 + hvx_vec_inverse_f32 avoids scalar extraction.
236    HVX_Vector sum_sf    = hvx_vec_reduce_sum_f32(Q6_Vsf_equals_Vqf32(sum_v));
237    HVX_Vector rsqrt_v   = hvx_vec_rsqrt_f32(sum_sf);              // 1/sqrt(sum)
238    HVX_Vector sqrt_v    = hvx_vec_inverse_f32(rsqrt_v);            // sqrt(sum)
239    HVX_Vector epsilon_v = hvx_vec_splat_f32(epsilon);
240    HVX_Vector denom_v   = Q6_Vsf_vmax_VsfVsf(sqrt_v, epsilon_v);  // fmax(sqrt(sum), epsilon)
241    HVX_Vector scale_v   = hvx_vec_inverse_f32(denom_v);            // 1/fmax(sqrt(sum), epsilon)
242
243    #pragma unroll(4)
244    for (int i = 0; i < nvec; i++) {
245        HVX_Vector v1 = v_src[i];
246        v_dst[i]      = Q6_Vsf_equals_Vqf32(Q6_Vqf32_vmpy_VsfVsf(v1, scale_v));
247    }
248
249    if (nloe > 0) {
250        HVX_VectorPred bmask = Q6_Q_vsetq_R(nloe * 4);
251        HVX_Vector v1 = Q6_V_vand_QV(bmask, v_src[nvec]);
252        HVX_Vector result = Q6_Vsf_equals_Vqf32(Q6_Vqf32_vmpy_VsfVsf(v1, scale_v));
253        hvx_vec_store_a(&v_dst[nvec], nloe * 4, result);
254    }
255}
256
257// F16 norm kernels: reduce and scale in f32 (via promote/narrow), matching the
258// precision-preserving pattern used by the flash-attn f16 kernels.
259
260static inline void hvx_fast_rms_norm_f16(const uint8_t * restrict src,
261                                         uint8_t * restrict dst,
262                                         const int num_elems,
263                                         float     epsilon) {
264
265    const HVX_Vector * restrict v_src = (HVX_Vector *) src;
266    HVX_Vector * restrict v_dst       = (HVX_Vector *) dst;
267
268    const int nvec = num_elems / VLEN_FP16;    // number of full f16 vectors
269    const int nloe = num_elems % VLEN_FP16;    // leftover elements
270
271    HVX_Vector sum_v = Q6_V_vsplat_R(0x00000000);
272    HVX_Vector epsilon_v = hvx_vec_splat_f32(epsilon);
273
274    #pragma unroll(4)
275    for (int i = 0; i < nvec; i++) {
276        HVX_VectorPair p = hvx_vec_f16_to_f32(v_src[i]);
277        HVX_Vector p0 = Q6_V_lo_W(p);
278        HVX_Vector p1 = Q6_V_hi_W(p);
279        sum_v = Q6_Vqf32_vadd_Vqf32Vqf32(sum_v, Q6_Vqf32_vmpy_VsfVsf(p0, p0));
280        sum_v = Q6_Vqf32_vadd_Vqf32Vqf32(sum_v, Q6_Vqf32_vmpy_VsfVsf(p1, p1));
281    }
282
283    if (nloe > 0) {
284        HVX_VectorPred bmask = Q6_Q_vsetq_R(nloe * SIZEOF_FP16);
285        HVX_Vector v1 = Q6_V_vand_QV(bmask, v_src[nvec]);
286        HVX_VectorPair p = hvx_vec_f16_to_f32(v1);
287        HVX_Vector p0 = Q6_V_lo_W(p);
288        HVX_Vector p1 = Q6_V_hi_W(p);
289        sum_v = Q6_Vqf32_vadd_Vqf32Vqf32(sum_v, Q6_Vqf32_vmpy_VsfVsf(p0, p0));
290        sum_v = Q6_Vqf32_vadd_Vqf32Vqf32(sum_v, Q6_Vqf32_vmpy_VsfVsf(p1, p1));
291    }
292
293    sum_v = hvx_vec_reduce_sum_f32(Q6_Vsf_equals_Vqf32(sum_v));
294
295    HVX_Vector t_v            = hvx_vec_splat_f32((float) num_elems);
296    HVX_Vector denom_v        = hvx_vec_inverse_f32(t_v);
297    HVX_Vector mean_v         = Q6_Vqf32_vmpy_VsfVsf(sum_v, denom_v);
298    HVX_Vector mean_epsilon_v = Q6_Vqf32_vadd_Vqf32Vsf(mean_v, epsilon_v);
299
300    HVX_Vector scale_v = hvx_vec_rsqrt_f32(Q6_Vsf_equals_Vqf32(mean_epsilon_v));
301
302    #pragma unroll(4)
303    for (int i = 0; i < nvec; i++) {
304        HVX_VectorPair p = hvx_vec_f16_to_f32(v_src[i]);
305        HVX_Vector r0 = Q6_Vsf_equals_Vqf32(Q6_Vqf32_vmpy_VsfVsf(Q6_V_lo_W(p), scale_v));
306        HVX_Vector r1 = Q6_Vsf_equals_Vqf32(Q6_Vqf32_vmpy_VsfVsf(Q6_V_hi_W(p), scale_v));
307        v_dst[i] = hvx_vec_f32_to_f16(r0, r1);
308    }
309
310    if (nloe > 0) {
311        HVX_VectorPred bmask = Q6_Q_vsetq_R(nloe * SIZEOF_FP16);
312        HVX_Vector v1 = Q6_V_vand_QV(bmask, v_src[nvec]);
313        HVX_VectorPair p = hvx_vec_f16_to_f32(v1);
314        HVX_Vector r0 = Q6_Vsf_equals_Vqf32(Q6_Vqf32_vmpy_VsfVsf(Q6_V_lo_W(p), scale_v));
315        HVX_Vector r1 = Q6_Vsf_equals_Vqf32(Q6_Vqf32_vmpy_VsfVsf(Q6_V_hi_W(p), scale_v));
316        HVX_Vector result = hvx_vec_f32_to_f16(r0, r1);
317        hvx_vec_store_a(&v_dst[nvec], nloe * SIZEOF_FP16, result);
318    }
319}
320
321static inline void hvx_fast_norm_f16(const uint8_t * restrict src,
322                                     uint8_t * restrict dst,
323                                     const int num_elems,
324                                     float     epsilon) {
325
326    const HVX_Vector * restrict v_src = (HVX_Vector *) src;
327    HVX_Vector * restrict v_dst       = (HVX_Vector *) dst;
328
329    const int nvec = num_elems / VLEN_FP16;
330    const int nloe = num_elems % VLEN_FP16;
331
332    HVX_Vector sum_sq_v = Q6_V_vsplat_R(0x00000000);
333    HVX_Vector sum_x_v  = Q6_V_vsplat_R(0x00000000);
334    HVX_Vector epsilon_v = hvx_vec_splat_f32(epsilon);
335
336    #pragma unroll(4)
337    for (int i = 0; i < nvec; i++) {
338        HVX_VectorPair p = hvx_vec_f16_to_f32(v_src[i]);
339        HVX_Vector p0 = Q6_V_lo_W(p);
340        HVX_Vector p1 = Q6_V_hi_W(p);
341        sum_sq_v = Q6_Vqf32_vadd_Vqf32Vqf32(sum_sq_v, Q6_Vqf32_vmpy_VsfVsf(p0, p0));
342        sum_sq_v = Q6_Vqf32_vadd_Vqf32Vqf32(sum_sq_v, Q6_Vqf32_vmpy_VsfVsf(p1, p1));
343        sum_x_v  = Q6_Vqf32_vadd_Vqf32Vqf32(sum_x_v,  Q6_Vqf32_vadd_VsfVsf(p0, Q6_V_vzero()));
344        sum_x_v  = Q6_Vqf32_vadd_Vqf32Vqf32(sum_x_v,  Q6_Vqf32_vadd_VsfVsf(p1, Q6_V_vzero()));
345    }
346
347    if (nloe > 0) {
348        HVX_VectorPred bmask = Q6_Q_vsetq_R(nloe * SIZEOF_FP16);
349        HVX_Vector v1 = Q6_V_vand_QV(bmask, v_src[nvec]);
350        HVX_VectorPair p = hvx_vec_f16_to_f32(v1);
351        HVX_Vector p0 = Q6_V_lo_W(p);
352        HVX_Vector p1 = Q6_V_hi_W(p);
353        sum_sq_v = Q6_Vqf32_vadd_Vqf32Vqf32(sum_sq_v, Q6_Vqf32_vmpy_VsfVsf(p0, p0));
354        sum_sq_v = Q6_Vqf32_vadd_Vqf32Vqf32(sum_sq_v, Q6_Vqf32_vmpy_VsfVsf(p1, p1));
355        sum_x_v  = Q6_Vqf32_vadd_Vqf32Vqf32(sum_x_v,  Q6_Vqf32_vadd_VsfVsf(p0, Q6_V_vzero()));
356        sum_x_v  = Q6_Vqf32_vadd_Vqf32Vqf32(sum_x_v,  Q6_Vqf32_vadd_VsfVsf(p1, Q6_V_vzero()));
357    }
358
359    sum_sq_v = hvx_vec_reduce_sum_f32(Q6_Vsf_equals_Vqf32(sum_sq_v));
360    sum_x_v  = hvx_vec_reduce_sum_f32(Q6_Vsf_equals_Vqf32(sum_x_v));
361
362    HVX_Vector t_v            = hvx_vec_splat_f32((float) num_elems);
363    HVX_Vector denom_v        = hvx_vec_inverse_f32(t_v);
364    HVX_Vector mean_sq_v      = Q6_Vqf32_vmpy_VsfVsf(sum_sq_v, denom_v);
365    HVX_Vector mean_x_v       = Q6_Vqf32_vmpy_VsfVsf(sum_x_v,  denom_v);
366    HVX_Vector mean_x_sq_v    = Q6_Vqf32_vmpy_VsfVsf(Q6_Vsf_equals_Vqf32(mean_x_v), Q6_Vsf_equals_Vqf32(mean_x_v));
367    HVX_Vector var_v          = Q6_Vqf32_vsub_Vqf32Vqf32(mean_sq_v, mean_x_sq_v);
368    HVX_Vector var_epsilon_v  = Q6_Vqf32_vadd_Vqf32Vsf(var_v, epsilon_v);
369
370    HVX_Vector scale_v  = hvx_vec_rsqrt_f32(Q6_Vsf_equals_Vqf32(var_epsilon_v));
371    HVX_Vector mean_x_b = hvx_vec_repl_f32(Q6_Vsf_equals_Vqf32(mean_x_v));
372
373    #pragma unroll(4)
374    for (int i = 0; i < nvec; i++) {
375        HVX_VectorPair p = hvx_vec_f16_to_f32(v_src[i]);
376        HVX_Vector d0 = Q6_Vqf32_vsub_VsfVsf(Q6_V_lo_W(p), mean_x_b);
377        HVX_Vector d1 = Q6_Vqf32_vsub_VsfVsf(Q6_V_hi_W(p), mean_x_b);
378        HVX_Vector r0 = Q6_Vsf_equals_Vqf32(Q6_Vqf32_vmpy_VsfVsf(Q6_Vsf_equals_Vqf32(d0), scale_v));
379        HVX_Vector r1 = Q6_Vsf_equals_Vqf32(Q6_Vqf32_vmpy_VsfVsf(Q6_Vsf_equals_Vqf32(d1), scale_v));
380        v_dst[i] = hvx_vec_f32_to_f16(r0, r1);
381    }
382
383    if (nloe > 0) {
384        HVX_VectorPred bmask = Q6_Q_vsetq_R(nloe * SIZEOF_FP16);
385        HVX_Vector v1 = Q6_V_vand_QV(bmask, v_src[nvec]);
386        HVX_VectorPair p = hvx_vec_f16_to_f32(v1);
387        HVX_Vector d0 = Q6_Vqf32_vsub_VsfVsf(Q6_V_lo_W(p), mean_x_b);
388        HVX_Vector d1 = Q6_Vqf32_vsub_VsfVsf(Q6_V_hi_W(p), mean_x_b);
389        HVX_Vector r0 = Q6_Vsf_equals_Vqf32(Q6_Vqf32_vmpy_VsfVsf(Q6_Vsf_equals_Vqf32(d0), scale_v));
390        HVX_Vector r1 = Q6_Vsf_equals_Vqf32(Q6_Vqf32_vmpy_VsfVsf(Q6_Vsf_equals_Vqf32(d1), scale_v));
391        HVX_Vector result = hvx_vec_f32_to_f16(r0, r1);
392        hvx_vec_store_a(&v_dst[nvec], nloe * SIZEOF_FP16, result);
393    }
394}
395
396static inline void hvx_fast_l2_norm_f16(const uint8_t * restrict src,
397                                        uint8_t * restrict dst,
398                                        const int num_elems,
399                                        float     epsilon) {
400
401    const HVX_Vector * restrict v_src = (HVX_Vector *) src;
402    HVX_Vector * restrict v_dst       = (HVX_Vector *) dst;
403
404    const int nvec = num_elems / VLEN_FP16;
405    const int nloe = num_elems % VLEN_FP16;
406
407    HVX_Vector sum_v = hvx_vec_splat_f32(0.0f);
408
409    #pragma unroll(4)
410    for (int i = 0; i < nvec; i++) {
411        HVX_VectorPair p = hvx_vec_f16_to_f32(v_src[i]);
412        HVX_Vector p0 = Q6_V_lo_W(p);
413        HVX_Vector p1 = Q6_V_hi_W(p);
414        sum_v = Q6_Vqf32_vadd_Vqf32Vqf32(sum_v, Q6_Vqf32_vmpy_VsfVsf(p0, p0));
415        sum_v = Q6_Vqf32_vadd_Vqf32Vqf32(sum_v, Q6_Vqf32_vmpy_VsfVsf(p1, p1));
416    }
417
418    if (nloe > 0) {
419        HVX_VectorPred bmask = Q6_Q_vsetq_R(nloe * SIZEOF_FP16);
420        HVX_Vector v1 = Q6_V_vand_QV(bmask, v_src[nvec]);
421        HVX_VectorPair p = hvx_vec_f16_to_f32(v1);
422        HVX_Vector p0 = Q6_V_lo_W(p);
423        HVX_Vector p1 = Q6_V_hi_W(p);
424        sum_v = Q6_Vqf32_vadd_Vqf32Vqf32(sum_v, Q6_Vqf32_vmpy_VsfVsf(p0, p0));
425        sum_v = Q6_Vqf32_vadd_Vqf32Vqf32(sum_v, Q6_Vqf32_vmpy_VsfVsf(p1, p1));
426    }
427
428    HVX_Vector sum_sf    = hvx_vec_reduce_sum_f32(Q6_Vsf_equals_Vqf32(sum_v));
429    HVX_Vector rsqrt_v   = hvx_vec_rsqrt_f32(sum_sf);
430    HVX_Vector sqrt_v    = hvx_vec_inverse_f32(rsqrt_v);
431    HVX_Vector epsilon_v = hvx_vec_splat_f32(epsilon);
432    HVX_Vector denom_v   = Q6_Vsf_vmax_VsfVsf(sqrt_v, epsilon_v);
433    HVX_Vector scale_v   = hvx_vec_inverse_f32(denom_v);
434
435    #pragma unroll(4)
436    for (int i = 0; i < nvec; i++) {
437        HVX_VectorPair p = hvx_vec_f16_to_f32(v_src[i]);
438        HVX_Vector r0 = Q6_Vsf_equals_Vqf32(Q6_Vqf32_vmpy_VsfVsf(Q6_V_lo_W(p), scale_v));
439        HVX_Vector r1 = Q6_Vsf_equals_Vqf32(Q6_Vqf32_vmpy_VsfVsf(Q6_V_hi_W(p), scale_v));
440        v_dst[i] = hvx_vec_f32_to_f16(r0, r1);
441    }
442
443    if (nloe > 0) {
444        HVX_VectorPred bmask = Q6_Q_vsetq_R(nloe * SIZEOF_FP16);
445        HVX_Vector v1 = Q6_V_vand_QV(bmask, v_src[nvec]);
446        HVX_VectorPair p = hvx_vec_f16_to_f32(v1);
447        HVX_Vector r0 = Q6_Vsf_equals_Vqf32(Q6_Vqf32_vmpy_VsfVsf(Q6_V_lo_W(p), scale_v));
448        HVX_Vector r1 = Q6_Vsf_equals_Vqf32(Q6_Vqf32_vmpy_VsfVsf(Q6_V_hi_W(p), scale_v));
449        HVX_Vector result = hvx_vec_f32_to_f16(r0, r1);
450        hvx_vec_store_a(&v_dst[nvec], nloe * SIZEOF_FP16, result);
451    }
452}
453
454#endif // HVX_NORM_H