char-fork/llama.cpp

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

Titaniumtownsycl: add a batched L2_NORM kernel (#28222)dbeb37548

char/aliceai-t5-moe
45.7 KiB1046 linesraw
1#include "norm.hpp"
2#include "ggml-sycl/common.hpp"
3#include "ggml-sycl/presets.hpp"
4
5static void norm_f32(const float* x, float* dst, const int ncols,
6    const int64_t src_stride_col, const int64_t src_stride_row, const int64_t src_stride_channel, const int64_t src_stride_sample,
7    const int64_t dst_stride_col, const int64_t dst_stride_row, const int64_t dst_stride_channel, const int64_t dst_stride_sample,
8    const float eps, const sycl::nd_item<3>& item_ct1, sycl::float2* s_sum, int block_size) {
9
10    const int nthreads = item_ct1.get_local_range(2);
11    const int sample  = item_ct1.get_group(0);
12    const int channel = item_ct1.get_group(1);
13    const int row     = item_ct1.get_group(2);
14
15    const int tid = item_ct1.get_local_id(2);
16    const int nwarps = nthreads / WARP_SIZE;
17
18    const auto src_offset = calculate_offset<3>({src_stride_sample, src_stride_channel, src_stride_row}, {sample, channel, row});
19    const auto dst_offset = calculate_offset<3>({dst_stride_sample, dst_stride_channel, dst_stride_row}, {sample, channel, row});
20
21    x += src_offset;
22    dst += dst_offset;
23
24    sycl::float2 mean_var = sycl::float2(0.f, 0.f);
25
26    for (int col = tid; col < ncols; col += block_size) {
27        const float xi = x[col * src_stride_col];
28        mean_var.x() += xi;
29        mean_var.y() += xi * xi;
30    }
31
32    // sum up partial sums
33    mean_var = warp_reduce_sum(mean_var, item_ct1);
34    if  (block_size > WARP_SIZE) {
35        const auto sub_group = item_ct1.get_sub_group();
36        const auto sg_id = sub_group.get_group_linear_id();
37        const auto wi_in_sg = sub_group.get_local_linear_id();
38        if (wi_in_sg == 0) {
39            s_sum[sg_id] = mean_var;
40        }
41        item_ct1.barrier(sycl::access::fence_space::local_space);
42        mean_var = 0.f;
43        const size_t nreduce = ceil_div(nwarps, WARP_SIZE);
44        for (size_t i = 0; i < nreduce; i += 1)
45        {
46            mean_var += s_sum[wi_in_sg + i * WARP_SIZE];
47        }
48        mean_var = warp_reduce_sum(mean_var, item_ct1);
49    }
50
51    const float mean = mean_var.x() / ncols;
52    const float var = mean_var.y() / ncols - mean * mean;
53    const float inv_std = sycl::rsqrt(var + eps);
54
55    for (int col = tid; col < ncols; col += block_size) {
56        dst[col * dst_stride_col] = (x[col * src_stride_col] - mean) * inv_std;
57    }
58}
59
60static void group_norm_f32(const float* x, float* dst, const int group_size, const int ne_elements, const float eps,
61    const sycl::nd_item<3>& item_ct1, float* s_sum, int block_size) {
62    int start = item_ct1.get_group(2) * group_size;
63    int end = start + group_size;
64    const int nthreads = item_ct1.get_local_range(2);
65    const int nwarps = nthreads / WARP_SIZE;
66    start += item_ct1.get_local_id(2);
67    size_t nreduce = nwarps / WARP_SIZE;
68
69    if (end >= ne_elements) {
70        end = ne_elements;
71    }
72
73    float tmp = 0.0f; // partial sum for thread in warp
74
75    for (int j = start; j < end; j += block_size) {
76        tmp += x[j];
77    }
78
79    tmp = warp_reduce_sum(tmp, item_ct1);
80    if (block_size > WARP_SIZE) {
81
82        int warp_id = item_ct1.get_local_id(2) / WARP_SIZE;
83        int lane_id = item_ct1.get_local_id(2) % WARP_SIZE;
84        if (lane_id == 0) {
85            s_sum[warp_id] = tmp;
86        }
87        /*
88        DPCT1118:1: SYCL group functions and algorithms must be encountered in
89        converged control flow. You may need to adjust the code.
90        */
91        /*
92        DPCT1065:54: Consider replacing sycl::nd_item::barrier() with
93        sycl::nd_item::barrier(sycl::access::fence_space::local_space) for
94        better performance if there is no access to global memory.
95        */
96        item_ct1.barrier();
97        tmp = 0.f;
98        for (size_t i = 0; i < nreduce; i += 1)
99        {
100            tmp += s_sum[lane_id + i * WARP_SIZE];
101        }
102        tmp = warp_reduce_sum(tmp, item_ct1);
103    }
104
105    float mean = tmp / group_size;
106    tmp = 0.0f;
107
108    for (int j = start; j < end; j += block_size) {
109        float xi = x[j] - mean;
110        dst[j] = xi;
111        tmp += xi * xi;
112    }
113
114    tmp = warp_reduce_sum(tmp, item_ct1);
115    if (block_size > WARP_SIZE) {
116
117        int warp_id = item_ct1.get_local_id(2) / WARP_SIZE;
118        int lane_id = item_ct1.get_local_id(2) % WARP_SIZE;
119        if (lane_id == 0) {
120            s_sum[warp_id] = tmp;
121        }
122        /*
123        DPCT1118:2: SYCL group functions and algorithms must be encountered in
124        converged control flow. You may need to adjust the code.
125        */
126        /*
127        DPCT1065:55: Consider replacing sycl::nd_item::barrier() with
128        sycl::nd_item::barrier(sycl::access::fence_space::local_space) for
129        better performance if there is no access to global memory.
130        */
131        item_ct1.barrier();
132        tmp = 0.f;
133        for (size_t i = 0; i < nreduce; i += 1)
134        {
135            tmp += s_sum[lane_id + i * WARP_SIZE];
136        }
137        tmp = warp_reduce_sum(tmp, item_ct1);
138    }
139
140    float variance = tmp / group_size;
141    float scale = sycl::rsqrt(variance + eps);
142    for (int j = start; j < end; j += block_size) {
143        dst[j] *= scale;
144    }
145}
146
147template <bool do_multiply = false, bool do_add = false>
148static void rms_norm_f32(const float* x, float* dst, const int ncols,
149    const int64_t src_stride_col, const int64_t src_stride_row, const int64_t src_stride_channel, const int64_t src_stride_sample,
150    const int64_t dst_stride_col, const int64_t dst_stride_row, const int64_t dst_stride_channel, const int64_t dst_stride_sample,
151    const float eps, const sycl::nd_item<3>& item_ct1, float* s_sum, int block_size,
152    const float* mul = nullptr, const int64_t mul_stride_row = 0, const int64_t mul_stride_channel = 0,
153    const int64_t mul_stride_sample = 0, const int mul_nrows = 0, const int mul_nchannels = 0, const int mul_nsamples = 0,
154    const float* add = nullptr, const int64_t add_stride_row = 0, const int64_t add_stride_channel = 0,
155    const int64_t add_stride_sample = 0, const int add_nrows = 0, const int add_nchannels = 0, const int add_nsamples = 0) {
156
157    static_assert(!do_add || do_multiply, "fusing add is not supported without multiplying");
158
159    const int sample  = item_ct1.get_group(0);
160    const int channel = item_ct1.get_group(1);
161    const int row     = item_ct1.get_group(2);
162
163    const int nthreads = item_ct1.get_local_range(2);
164
165    const int tid = item_ct1.get_local_id(2);
166    const int nwarps = nthreads / WARP_SIZE;
167
168    const auto src_offset = calculate_offset<3>({src_stride_sample, src_stride_channel, src_stride_row}, {sample, channel, row});
169    const auto dst_offset = calculate_offset<3>({dst_stride_sample, dst_stride_channel, dst_stride_row}, {sample, channel, row});
170
171    x   += src_offset;
172    dst += dst_offset;
173
174    if constexpr (do_multiply) {
175        const int mul_row     = row     % mul_nrows;
176        const int mul_channel = channel % mul_nchannels;
177        const int mul_sample  = sample  % mul_nsamples;
178        mul += mul_sample * mul_stride_sample + mul_channel * mul_stride_channel + mul_row * mul_stride_row;
179    }
180
181    if constexpr (do_add) {
182        const int add_row     = row     % add_nrows;
183        const int add_channel = channel % add_nchannels;
184        const int add_sample  = sample  % add_nsamples;
185        add += add_sample * add_stride_sample + add_channel * add_stride_channel + add_row * add_stride_row;
186    }
187
188    float tmp = 0.0f; // partial sum for thread in warp
189
190    for (int col = tid; col < ncols; col += block_size) {
191        const float xi = x[col * src_stride_col];
192        tmp += xi * xi;
193    }
194
195    // sum up partial sums
196    tmp = warp_reduce_sum(tmp, item_ct1);
197    if (block_size > WARP_SIZE) {
198        const auto sub_group = item_ct1.get_sub_group();
199        const auto sg_id = sub_group.get_group_linear_id();
200        const auto wi_in_sg = sub_group.get_local_linear_id();
201        if (wi_in_sg == 0) {
202            s_sum[sg_id] = tmp;
203        }
204
205        item_ct1.barrier(sycl::access::fence_space::local_space);
206        const size_t nreduce = ceil_div(nwarps, WARP_SIZE);
207        tmp = 0.f;
208        for (size_t i = 0; i < nreduce; i += 1)
209        {
210            tmp += s_sum[wi_in_sg + i * WARP_SIZE];
211        }
212        tmp = warp_reduce_sum(tmp, item_ct1);
213    }
214
215    const float mean = tmp / ncols;
216    const float scale = sycl::rsqrt(mean + eps);
217
218    for (int col = tid; col < ncols; col += block_size) {
219        if constexpr (do_multiply && do_add) {
220            dst[col * dst_stride_col] = scale * x[col * src_stride_col] * mul[col] + add[col];
221        } else if constexpr (do_multiply) {
222            dst[col * dst_stride_col] = scale * x[col * src_stride_col] * mul[col];
223        } else {
224            dst[col * dst_stride_col] = scale * x[col * src_stride_col];
225        }
226    }
227}
228
229template<int warp_size>
230static void l2_norm_f32(const float * x, float * dst, const int ncols,
231    const int64_t src_stride_col, const int64_t src_stride_row, const int64_t src_stride_channel,
232    const int64_t src_stride_sample, const int64_t dst_stride_col, const int64_t dst_stride_row,
233    const int64_t dst_stride_channel, const int64_t dst_stride_sample, const float eps,
234    const sycl::nd_item<3>& item_ct1, float* s_sum, const int block_size) {
235
236    const int row     = item_ct1.get_group(2);
237    const int channel = item_ct1.get_group(1);
238    const int sample  = item_ct1.get_group(0);
239    const int tid     = item_ct1.get_local_id(2);
240
241    x   += sample*src_stride_sample + channel*src_stride_channel + row*src_stride_row;
242    dst += sample*dst_stride_sample + channel*dst_stride_channel + row*dst_stride_row;
243
244    float tmp = 0.0f; // partial sum for thread in warp
245
246    for (int col = tid; col < ncols; col += block_size) {
247        const float xi = x[col * src_stride_col];
248        tmp += xi * xi;
249    }
250
251    tmp = block_reduce<block_reduce_method::SUM, warp_size>(tmp, s_sum, block_size);
252    const float scale = sycl::rsqrt(sycl::fmax(tmp, eps * eps));
253
254    for (int col = tid; col < ncols; col += block_size) {
255        dst[col * dst_stride_col] = scale * x[col * src_stride_col];
256    }
257}
258
259static void norm_f32_sycl(const float * x, float * dst, const int ncols, const int nrows, const int nchannels, const int nsamples,
260    const int64_t src_stride_col, const int64_t src_stride_row, const int64_t src_stride_channel, const int64_t src_stride_sample,
261    const int64_t dst_stride_col, const int64_t dst_stride_row, const int64_t dst_stride_channel, const int64_t dst_stride_sample,
262        const float eps, queue_ptr stream, int device) {
263
264    const sycl::range<3> global_dims(nsamples, nchannels, nrows);
265    if (ncols < 1024) {
266        const sycl::range<3> block_dims(1, 1, WARP_SIZE);
267        stream->submit([&](sycl::handler& cgh) {
268            cgh.parallel_for(
269                sycl::nd_range<3>(global_dims * block_dims, block_dims),
270                [=](sycl::nd_item<3> item_ct1)
271                [[sycl::reqd_sub_group_size(WARP_SIZE)]] {
272                    norm_f32(x, dst, ncols,
273                        src_stride_col, src_stride_row, src_stride_channel, src_stride_sample,
274                        dst_stride_col, dst_stride_row, dst_stride_channel, dst_stride_sample,
275                        eps, item_ct1, nullptr, WARP_SIZE);
276                });
277            });
278    }
279    else {
280        const int work_group_size = ggml_sycl_info().max_work_group_sizes[device];
281        assert(work_group_size % (WARP_SIZE * WARP_SIZE) == 0);
282        const sycl::range<3> block_dims(1, 1, work_group_size);
283        /*
284        DPCT1049:17: The work-group size passed to the SYCL kernel may exceed
285        the limit. To get the device limit, query
286        info::device::max_work_group_size. Adjust the work-group size if needed.
287        */
288        stream->submit([&](sycl::handler& cgh) {
289            sycl::local_accessor<sycl::float2, 1> s_sum_acc_ct1(
290                            sycl::range<1>(work_group_size / WARP_SIZE), cgh);
291            cgh.parallel_for(
292                sycl::nd_range<3>(global_dims * block_dims, block_dims),
293                [=](sycl::nd_item<3> item_ct1)
294                [[sycl::reqd_sub_group_size(WARP_SIZE)]] {
295                    norm_f32(x, dst, ncols,
296                        src_stride_col, src_stride_row, src_stride_channel, src_stride_sample,
297                        dst_stride_col, dst_stride_row, dst_stride_channel, dst_stride_sample,
298                        eps, item_ct1, get_pointer(s_sum_acc_ct1), work_group_size);
299                });
300            });
301    }
302}
303
304static void group_norm_f32_sycl(const float* x, float* dst,
305    const int num_groups, const float eps, const int group_size,
306    const int ne_elements, queue_ptr stream, int device) {
307    if (group_size < 1024) {
308        const sycl::range<3> block_dims(1, 1, WARP_SIZE);
309        stream->submit([&](sycl::handler& cgh) {
310            const float eps_ct4 = eps;
311            cgh.parallel_for(
312                sycl::nd_range<3>(sycl::range<3>(1, 1, num_groups) * block_dims,
313                    block_dims),
314                [=](sycl::nd_item<3> item_ct1)
315                [[sycl::reqd_sub_group_size(WARP_SIZE)]] {
316                    group_norm_f32(
317                        x, dst, group_size, ne_elements, eps_ct4, item_ct1,
318                        nullptr, WARP_SIZE);
319                });
320            });
321    }
322    else {
323        const int work_group_size = ggml_sycl_info().max_work_group_sizes[device];
324        assert(work_group_size % (WARP_SIZE * WARP_SIZE) == 0);
325        const sycl::range<3> block_dims(1, 1, work_group_size);
326        /*
327        DPCT1049:18: The work-group size passed to the SYCL kernel may exceed
328        the limit. To get the device limit, query
329        info::device::max_work_group_size. Adjust the work-group size if needed.
330        */
331
332        stream->submit([&](sycl::handler& cgh) {
333            sycl::local_accessor<float, 1> s_sum_acc_ct1(sycl::range<1>(work_group_size / WARP_SIZE),
334                cgh);
335
336            const float eps_ct4 = eps;
337
338            cgh.parallel_for(
339                sycl::nd_range<3>(sycl::range<3>(1, 1, num_groups) * block_dims,
340                    block_dims),
341                [=](sycl::nd_item<3> item_ct1)
342                [[sycl::reqd_sub_group_size(WARP_SIZE)]] {
343                    group_norm_f32(x, dst, group_size, ne_elements,
344                        eps_ct4, item_ct1,
345                        get_pointer(s_sum_acc_ct1), work_group_size);
346                });
347            });
348    }
349}
350
351static void rms_norm_f32_sycl(const float* x, float* dst, const int ncols, const int nrows, const int nchannels, const int nsamples,
352    const int64_t src_stride_col, const int64_t src_stride_row, const int64_t src_stride_channel, const int64_t src_stride_sample,
353    const int64_t dst_stride_col, const int64_t dst_stride_row, const int64_t dst_stride_channel, const int64_t dst_stride_sample,
354    const float eps, queue_ptr stream, int device) {
355    // printf("%s ncols=%d, nrows=%d, WARP_SIZE=%d\n", __func__, ncols, nrows, WARP_SIZE);
356
357    const sycl::range<3> global_dims(nsamples, nchannels, nrows);
358    if (ncols < 1024) {
359        const sycl::range<3> block_dims(1, 1, WARP_SIZE);
360        stream->submit([&](sycl::handler& cgh) {
361            cgh.parallel_for(
362                sycl::nd_range<3>(global_dims * block_dims, block_dims),
363                [=](sycl::nd_item<3> item_ct1)
364                [[sycl::reqd_sub_group_size(WARP_SIZE)]] {
365                    rms_norm_f32(x, dst, ncols,
366                        src_stride_col, src_stride_row, src_stride_channel, src_stride_sample,
367                        dst_stride_col, dst_stride_row, dst_stride_channel, dst_stride_sample,
368                        eps, item_ct1, nullptr, WARP_SIZE);
369                });
370            });
371    }
372    else {
373        const int work_group_size = ggml_sycl_info().max_work_group_sizes[device];
374        assert(work_group_size % (WARP_SIZE * WARP_SIZE) == 0);
375        const sycl::range<3> block_dims(1, 1, work_group_size);
376        /*
377        DPCT1049:19: The work-group size passed to the SYCL kernel may exceed
378        the limit. To get the device limit, query
379        info::device::max_work_group_size. Adjust the work-group size if needed.
380        */
381        stream->submit([&](sycl::handler& cgh) {
382            sycl::local_accessor<float, 1> s_sum_acc_ct1(sycl::range<1>(work_group_size / WARP_SIZE),
383                cgh);
384            cgh.parallel_for(
385                sycl::nd_range<3>(global_dims * block_dims, block_dims),
386                [=](sycl::nd_item<3> item_ct1)
387                [[sycl::reqd_sub_group_size(WARP_SIZE)]] {
388                    rms_norm_f32(x, dst, ncols,
389                        src_stride_col, src_stride_row, src_stride_channel, src_stride_sample,
390                        dst_stride_col, dst_stride_row, dst_stride_channel, dst_stride_sample,
391                        eps, item_ct1, get_pointer(s_sum_acc_ct1), work_group_size);
392                });
393            });
394    }
395}
396
397static void rms_norm_mul_f32_sycl(const float* x, const float* mul, float* dst, const int ncols, const int nrows,
398        const int nchannels, const int nsamples,
399        const int64_t src_stride_col, const int64_t src_stride_row, const int64_t src_stride_channel, const int64_t src_stride_sample,
400        const int64_t dst_stride_col, const int64_t dst_stride_row, const int64_t dst_stride_channel, const int64_t dst_stride_sample,
401        const int64_t mul_stride_row, const int64_t mul_stride_channel, const int64_t mul_stride_sample,
402        const int mul_nrows, const int mul_nchannels, const int mul_nsamples,
403        const float eps, queue_ptr stream, int device) {
404    const sycl::range<3> global_dims(nsamples, nchannels, nrows);
405    if (ncols < 1024) {
406        const sycl::range<3> block_dims(1, 1, WARP_SIZE);
407        stream->submit([&](sycl::handler& cgh) {
408            cgh.parallel_for(
409                sycl::nd_range<3>(global_dims * block_dims, block_dims),
410                [=](sycl::nd_item<3> item_ct1)
411                [[sycl::reqd_sub_group_size(WARP_SIZE)]] {
412                    rms_norm_f32<true>(x, dst, ncols,
413                        src_stride_col, src_stride_row, src_stride_channel, src_stride_sample,
414                        dst_stride_col, dst_stride_row, dst_stride_channel, dst_stride_sample,
415                        eps, item_ct1, nullptr, WARP_SIZE,
416                        mul, mul_stride_row, mul_stride_channel, mul_stride_sample, mul_nrows, mul_nchannels, mul_nsamples);
417                });
418            });
419    }
420    else {
421        const int work_group_size = ggml_sycl_info().max_work_group_sizes[device];
422        assert(work_group_size % (WARP_SIZE * WARP_SIZE) == 0);
423        const sycl::range<3> block_dims(1, 1, work_group_size);
424        stream->submit([&](sycl::handler& cgh) {
425            sycl::local_accessor<float, 1> s_sum_acc_ct1(sycl::range<1>(work_group_size / WARP_SIZE), cgh);
426            cgh.parallel_for(
427                sycl::nd_range<3>(global_dims * block_dims, block_dims),
428                [=](sycl::nd_item<3> item_ct1)
429                [[sycl::reqd_sub_group_size(WARP_SIZE)]] {
430                    rms_norm_f32<true>(x, dst, ncols,
431                        src_stride_col, src_stride_row, src_stride_channel, src_stride_sample,
432                        dst_stride_col, dst_stride_row, dst_stride_channel, dst_stride_sample,
433                        eps, item_ct1, get_pointer(s_sum_acc_ct1), work_group_size,
434                        mul, mul_stride_row, mul_stride_channel, mul_stride_sample, mul_nrows, mul_nchannels, mul_nsamples);
435                });
436            });
437    }
438}
439
440static void rms_norm_mul_add_f32_sycl(const float* x, const float* mul, const float* add, float* dst,
441        const int ncols, const int nrows, const int nchannels, const int nsamples,
442        const int64_t src_stride_col, const int64_t src_stride_row, const int64_t src_stride_channel, const int64_t src_stride_sample,
443        const int64_t dst_stride_col, const int64_t dst_stride_row, const int64_t dst_stride_channel, const int64_t dst_stride_sample,
444        const int64_t mul_stride_row, const int64_t mul_stride_channel, const int64_t mul_stride_sample,
445        const int mul_nrows, const int mul_nchannels, const int mul_nsamples,
446        const int64_t add_stride_row, const int64_t add_stride_channel, const int64_t add_stride_sample,
447        const int add_nrows, const int add_nchannels, const int add_nsamples,
448        const float eps, queue_ptr stream, int device) {
449    const sycl::range<3> global_dims(nsamples, nchannels, nrows);
450    if (ncols < 1024) {
451        const sycl::range<3> block_dims(1, 1, WARP_SIZE);
452        stream->submit([&](sycl::handler& cgh) {
453            cgh.parallel_for(
454                sycl::nd_range<3>(global_dims * block_dims, block_dims),
455                [=](sycl::nd_item<3> item_ct1)
456                [[sycl::reqd_sub_group_size(WARP_SIZE)]] {
457                    rms_norm_f32<true, true>(x, dst, ncols,
458                        src_stride_col, src_stride_row, src_stride_channel, src_stride_sample,
459                        dst_stride_col, dst_stride_row, dst_stride_channel, dst_stride_sample,
460                        eps, item_ct1, nullptr, WARP_SIZE,
461                        mul, mul_stride_row, mul_stride_channel, mul_stride_sample, mul_nrows, mul_nchannels, mul_nsamples,
462                        add, add_stride_row, add_stride_channel, add_stride_sample, add_nrows, add_nchannels, add_nsamples);
463                });
464            });
465    }
466    else {
467        const int work_group_size = ggml_sycl_info().max_work_group_sizes[device];
468        assert(work_group_size % (WARP_SIZE * WARP_SIZE) == 0);
469        const sycl::range<3> block_dims(1, 1, work_group_size);
470        stream->submit([&](sycl::handler& cgh) {
471            sycl::local_accessor<float, 1> s_sum_acc_ct1(sycl::range<1>(work_group_size / WARP_SIZE), cgh);
472            cgh.parallel_for(
473                sycl::nd_range<3>(global_dims * block_dims, block_dims),
474                [=](sycl::nd_item<3> item_ct1)
475                [[sycl::reqd_sub_group_size(WARP_SIZE)]] {
476                    rms_norm_f32<true, true>(x, dst, ncols,
477                        src_stride_col, src_stride_row, src_stride_channel, src_stride_sample,
478                        dst_stride_col, dst_stride_row, dst_stride_channel, dst_stride_sample,
479                        eps, item_ct1, get_pointer(s_sum_acc_ct1), work_group_size,
480                        mul, mul_stride_row, mul_stride_channel, mul_stride_sample, mul_nrows, mul_nchannels, mul_nsamples,
481                        add, add_stride_row, add_stride_channel, add_stride_sample, add_nrows, add_nchannels, add_nsamples);
482                });
483            });
484    }
485}
486
487template<int warp_size>
488static void l2_norm_f32_sycl(const float *   x,
489                             float *         dst,
490                             const int       ncols,
491                             const int       nrows,
492                             const int       nchannels,
493                             const int       nsamples,
494                             const int64_t   src_stride_col,
495                             const int64_t   src_stride_row,
496                             const int64_t   src_stride_channel,
497                             const int64_t   src_stride_sample,
498                             const int64_t   dst_stride_col,
499                             const int64_t   dst_stride_row,
500                             const int64_t   dst_stride_channel,
501                             const int64_t   dst_stride_sample,
502                             const float     eps,
503                             queue_ptr       stream,
504                             int             device) {
505    const dpct::dim3 blocks_num(nrows, nchannels, nsamples);
506
507    if (ncols < 1024) {
508        const dpct::dim3 block_dims(warp_size, 1, 1);
509        stream->submit([&](sycl::handler& cgh) {
510            cgh.parallel_for(
511                sycl::nd_range<3>(blocks_num * block_dims,
512                    block_dims),
513                [=](sycl::nd_item<3> item_ct1)
514                [[sycl::reqd_sub_group_size(warp_size)]] {
515                    l2_norm_f32<warp_size>(x, dst, ncols,
516                        src_stride_col, src_stride_row, src_stride_channel, src_stride_sample,
517                        dst_stride_col, dst_stride_row, dst_stride_channel, dst_stride_sample,
518                        eps, item_ct1,
519                        nullptr, warp_size);
520                });
521            });
522    }
523    else {
524        const int work_group_size = ggml_sycl_info().max_work_group_sizes[device];
525        assert(work_group_size % (warp_size * warp_size) == 0);
526        const sycl::range<3> block_dims(1, 1, work_group_size);
527        int lsm_size =  block_dims[2] > warp_size ? work_group_size / warp_size * sizeof(float): 0;
528        stream->submit([&](sycl::handler& cgh) {
529            sycl::local_accessor<float, 1> s_sum_acc_ct1(sycl::range<1>(lsm_size),
530                cgh);
531
532            cgh.parallel_for(
533                sycl::nd_range<3>(blocks_num * block_dims,
534                    block_dims),
535                [=](sycl::nd_item<3> item_ct1)
536                [[sycl::reqd_sub_group_size(warp_size)]] {
537                    l2_norm_f32<warp_size>(x, dst, ncols,
538                        src_stride_col, src_stride_row, src_stride_channel, src_stride_sample,
539                        dst_stride_col, dst_stride_row, dst_stride_channel, dst_stride_sample,
540                        eps, item_ct1, get_pointer(s_sum_acc_ct1), work_group_size);
541                });
542            });
543    }
544}
545
546// Batched L2 norm: N independent same-shape F32 tensors in one launch; the tensor
547// index is folded into grid dim0 and each row's reduction is identical to the
548// single-tensor kernel, so the result is bit-exact.
549struct l2_batch_ptrs {
550    const float * src[GGML_SYCL_L2_BATCH_MAX];
551    float *       dst[GGML_SYCL_L2_BATCH_MAX];
552};
553
554// One stride set shared by the whole batch: the caller only groups tensors whose nb[]
555// all match, so per-tensor state stays two pointers.
556struct l2_batch_strides {
557    int     ne1, ne2;
558    int64_t ss0, ss1, ss2, ss3;
559    int64_t ds0, ds1, ds2, ds3;
560};
561
562template <int warp_size>
563static void l2_norm_f32_batch(l2_batch_ptrs p, l2_batch_strides st, const int ncols, const float eps,
564                              const sycl::nd_item<3> & item_ct1) {
565    const int t   = item_ct1.get_group(0);  // tensor index
566    const int r   = item_ct1.get_group(2);  // flattened row over ne1*ne2*ne3
567    const int tid = item_ct1.get_local_id(2);
568
569    const int i1 = r % st.ne1;
570    const int i2 = (r / st.ne1) % st.ne2;
571    const int i3 = r / (st.ne1 * st.ne2);
572
573    const float * x   = p.src[t] + i3 * st.ss3 + i2 * st.ss2 + i1 * st.ss1;
574    float *       dst = p.dst[t] + i3 * st.ds3 + i2 * st.ds2 + i1 * st.ds1;
575
576    float tmp = 0.0f;
577    for (int col = tid; col < ncols; col += warp_size) {
578        const float xi = x[col * st.ss0];
579        tmp += xi * xi;
580    }
581    tmp = block_reduce<block_reduce_method::SUM, warp_size>(tmp, (float *) nullptr, warp_size);
582    const float scale = sycl::rsqrt(sycl::fmax(tmp, eps * eps));
583    for (int col = tid; col < ncols; col += warp_size) {
584        dst[col * st.ds0] = scale * x[col * st.ss0];
585    }
586}
587
588template <int warp_size>
589static void l2_norm_f32_batch_sycl(l2_batch_ptrs p, l2_batch_strides st, const int n_tensors,
590                                   const int ncols, const int nrows_total, const float eps,
591                                   queue_ptr stream) {
592    const dpct::dim3 blocks_num(nrows_total, 1, n_tensors);
593    const dpct::dim3 block_dims(warp_size, 1, 1);
594    stream->submit([&](sycl::handler & cgh) {
595        cgh.parallel_for(sycl::nd_range<3>(blocks_num * block_dims, block_dims),
596            [=](sycl::nd_item<3> item_ct1) [[sycl::reqd_sub_group_size(warp_size)]] {
597                l2_norm_f32_batch<warp_size>(p, st, ncols, eps, item_ct1);
598            });
599    });
600}
601
602void ggml_sycl_op_norm(ggml_backend_sycl_context& ctx, ggml_tensor* dst) {
603    const ggml_tensor * src0 = dst->src[0];
604
605    GGML_ASSERT(dst->src[0]->type == GGML_TYPE_F32);
606    GGML_ASSERT(dst->type == GGML_TYPE_F32);
607
608    GGML_TENSOR_UNARY_OP_LOCALS
609    dpct::queue_ptr main_stream = ctx.stream();
610    SYCL_CHECK(ggml_sycl_set_device(ctx.device));
611    const float * src0_dd = static_cast<const float *>(dst->src[0]->data);
612    float *       dst_dd  = static_cast<float *>(dst->data);
613
614    float eps;
615    memcpy(&eps, dst->op_params, sizeof(float));
616    GGML_ASSERT(eps >= 0.0f);
617    const size_t ts0 = ggml_type_size(src0->type);
618    const size_t tdst = ggml_type_size(dst->type);
619    GGML_ASSERT(nb00 % ts0 == 0 && nb01 % ts0 == 0 && nb02 % ts0 == 0 && nb03 % ts0 == 0);
620    GGML_ASSERT(nb0 % tdst == 0 && nb1 % tdst == 0 && nb2 % tdst == 0 && nb3 % tdst == 0);
621    const int64_t ss0 = nb00 / ts0;
622    const int64_t ss1 = nb01 / ts0;
623    const int64_t ss2 = nb02 / ts0;
624    const int64_t ss3 = nb03 / ts0;
625    const int64_t ds0 = nb0 / tdst;
626    const int64_t ds1 = nb1 / tdst;
627    const int64_t ds2 = nb2 / tdst;
628    const int64_t ds3 = nb3 / tdst;
629
630    norm_f32_sycl(src0_dd, dst_dd, ne00, ne01, ne02, ne03,
631        ss0, ss1, ss2, ss3, ds0, ds1, ds2, ds3, eps, main_stream, ctx.device);
632}
633
634void ggml_sycl_op_group_norm(ggml_backend_sycl_context& ctx, ggml_tensor* dst) {
635
636    GGML_ASSERT(dst->src[0]->type == GGML_TYPE_F32);
637    GGML_ASSERT(dst->type == GGML_TYPE_F32);
638
639    int num_groups = dst->op_params[0];
640    dpct::queue_ptr main_stream = ctx.stream();
641    SYCL_CHECK(ggml_sycl_set_device(ctx.device));
642
643    const float * src0_dd = static_cast<const float *>(dst->src[0]->data);
644    float *       dst_dd  = static_cast<float *>(dst->data);
645
646    float eps;
647    memcpy(&eps, dst->op_params + 1, sizeof(float));
648
649    int group_size = dst->src[0]->ne[0] * dst->src[0]->ne[1] * ((dst->src[0]->ne[2] + num_groups - 1) / num_groups);
650    group_norm_f32_sycl(src0_dd, dst_dd, num_groups, eps, group_size, dst->src[0]->ne[0] * dst->src[0]->ne[1] * dst->src[0]->ne[2], main_stream, ctx.device);
651}
652
653void ggml_sycl_op_rms_norm(ggml_backend_sycl_context & ctx, ggml_tensor * dst) {
654
655    const ggml_tensor * src0 = dst->src[0];
656    GGML_ASSERT(dst->src[0]->type == GGML_TYPE_F32);
657    GGML_ASSERT(dst->type == GGML_TYPE_F32);
658
659    dpct::queue_ptr main_stream = ctx.stream();
660    SYCL_CHECK(ggml_sycl_set_device(ctx.device));
661
662    const float * src0_dd = static_cast<const float *>(dst->src[0]->data);
663    float *       dst_dd  = static_cast<float *>(dst->data);
664
665    float eps;
666    memcpy(&eps, dst->op_params, sizeof(float));
667
668    GGML_TENSOR_UNARY_OP_LOCALS
669    const size_t ts0 = ggml_type_size(src0->type);
670    const size_t tdst = ggml_type_size(dst->type);
671    GGML_ASSERT(nb00 % ts0 == 0 && nb01 % ts0 == 0 && nb02 % ts0 == 0 && nb03 % ts0 == 0);
672    GGML_ASSERT(nb0 % tdst == 0 && nb1 % tdst == 0 && nb2 % tdst == 0 && nb3 % tdst == 0);
673    const int64_t ss0 = nb00 / ts0;
674    const int64_t ss1 = nb01 / ts0;
675    const int64_t ss2 = nb02 / ts0;
676    const int64_t ss3 = nb03 / ts0;
677    const int64_t ds0 = nb0 / tdst;
678    const int64_t ds1 = nb1 / tdst;
679    const int64_t ds2 = nb2 / tdst;
680    const int64_t ds3 = nb3 / tdst;
681    rms_norm_f32_sycl(src0_dd, dst_dd, ne00, ne01, ne02, ne03,
682        ss0, ss1, ss2, ss3, ds0, ds1, ds2, ds3, eps, main_stream, ctx.device);
683}
684
685void ggml_sycl_op_rms_norm_fused(ggml_backend_sycl_context & ctx, ggml_tensor * dst, ggml_tensor * mul_tensor) {
686    const ggml_tensor * rms_norm_src = dst->src[0];
687    float eps = 0.0f;
688    memcpy(&eps, dst->op_params, sizeof(float));
689
690    const float *       src0_dd = static_cast<const float *>(rms_norm_src->data);
691    const float *       mul_dd  = nullptr;
692    const ggml_tensor * mul_src = nullptr;
693    if (mul_tensor->src[0] == dst) {
694        mul_dd  = static_cast<const float *>(mul_tensor->src[1]->data);
695        mul_src = mul_tensor->src[1];
696    } else if (mul_tensor->src[1] == dst) {
697        mul_dd  = static_cast<const float *>(mul_tensor->src[0]->data);
698        mul_src = mul_tensor->src[0];
699    } else {
700        GGML_ASSERT(false);
701    }
702    float * dst_dd = static_cast<float *>(mul_tensor->data);
703
704    dpct::queue_ptr main_stream = ctx.stream();
705    SYCL_CHECK(ggml_sycl_set_device(ctx.device));
706
707    GGML_ASSERT(rms_norm_src->type == GGML_TYPE_F32);
708    GGML_ASSERT(dst->type == GGML_TYPE_F32);
709    GGML_ASSERT(mul_tensor->type == GGML_TYPE_F32);
710    GGML_ASSERT(eps >= 0.0f);
711
712    const int64_t ne00 = rms_norm_src->ne[0];
713    const int64_t ne01 = rms_norm_src->ne[1];
714    const int64_t ne02 = rms_norm_src->ne[2];
715    const int64_t ne03 = rms_norm_src->ne[3];
716
717    const size_t ts0 = ggml_type_size(rms_norm_src->type);
718    GGML_ASSERT(rms_norm_src->nb[0] == ts0);
719    const int64_t s00 = rms_norm_src->nb[0] / ts0;
720    const int64_t s01 = rms_norm_src->nb[1] / ts0;
721    const int64_t s02 = rms_norm_src->nb[2] / ts0;
722    const int64_t s03 = rms_norm_src->nb[3] / ts0;
723
724    const size_t tdst = ggml_type_size(mul_tensor->type);
725    GGML_ASSERT(mul_tensor->nb[0] == tdst);
726    const int64_t d00 = mul_tensor->nb[0] / tdst;
727    const int64_t d01 = mul_tensor->nb[1] / tdst;
728    const int64_t d02 = mul_tensor->nb[2] / tdst;
729    const int64_t d03 = mul_tensor->nb[3] / tdst;
730
731    const size_t ts_mul = ggml_type_size(mul_src->type);
732    GGML_ASSERT(mul_src->nb[0] == ts_mul);
733    const int64_t mul_s01 = mul_src->nb[1] / ts_mul;
734    const int64_t mul_s02 = mul_src->nb[2] / ts_mul;
735    const int64_t mul_s03 = mul_src->nb[3] / ts_mul;
736    const int mul_nrows     = mul_src->ne[1];
737    const int mul_nchannels = mul_src->ne[2];
738    const int mul_nsamples  = mul_src->ne[3];
739
740    rms_norm_mul_f32_sycl(src0_dd, mul_dd, dst_dd, ne00, ne01, ne02, ne03,
741        s00, s01, s02, s03, d00, d01, d02, d03,
742        mul_s01, mul_s02, mul_s03, mul_nrows, mul_nchannels, mul_nsamples, eps, main_stream, ctx.device);
743}
744
745void ggml_sycl_op_rms_norm_fused_add(ggml_backend_sycl_context & ctx, ggml_tensor * dst,
746                                     ggml_tensor * mul_tensor, ggml_tensor * add_tensor) {
747    const ggml_tensor * rms_norm_src = dst->src[0];
748    float eps = 0.0f;
749    memcpy(&eps, dst->op_params, sizeof(float));
750
751    const float *       src0_dd = static_cast<const float *>(rms_norm_src->data);
752    const float *       mul_dd  = nullptr;
753    const ggml_tensor * mul_src = nullptr;
754    if (mul_tensor->src[0] == dst) {
755        mul_dd  = static_cast<const float *>(mul_tensor->src[1]->data);
756        mul_src = mul_tensor->src[1];
757    } else if (mul_tensor->src[1] == dst) {
758        mul_dd  = static_cast<const float *>(mul_tensor->src[0]->data);
759        mul_src = mul_tensor->src[0];
760    } else {
761        GGML_ASSERT(false);
762    }
763
764    const float *       add_dd  = nullptr;
765    const ggml_tensor * add_src = nullptr;
766    if (add_tensor->src[0] == mul_tensor) {
767        add_dd  = static_cast<const float *>(add_tensor->src[1]->data);
768        add_src = add_tensor->src[1];
769    } else if (add_tensor->src[1] == mul_tensor) {
770        add_dd  = static_cast<const float *>(add_tensor->src[0]->data);
771        add_src = add_tensor->src[0];
772    } else {
773        GGML_ASSERT(false);
774    }
775
776    float * dst_dd = static_cast<float *>(add_tensor->data);
777
778    dpct::queue_ptr main_stream = ctx.stream();
779    SYCL_CHECK(ggml_sycl_set_device(ctx.device));
780
781    GGML_ASSERT(rms_norm_src->type == GGML_TYPE_F32);
782    GGML_ASSERT(dst->type == GGML_TYPE_F32);
783    GGML_ASSERT(mul_tensor->type == GGML_TYPE_F32);
784    GGML_ASSERT(add_tensor->type == GGML_TYPE_F32);
785    GGML_ASSERT(eps >= 0.0f);
786
787    const int64_t ne00 = rms_norm_src->ne[0];
788    const int64_t ne01 = rms_norm_src->ne[1];
789    const int64_t ne02 = rms_norm_src->ne[2];
790    const int64_t ne03 = rms_norm_src->ne[3];
791
792    const size_t ts0 = ggml_type_size(rms_norm_src->type);
793    GGML_ASSERT(rms_norm_src->nb[0] == ts0);
794    const int64_t s00 = rms_norm_src->nb[0] / ts0;
795    const int64_t s01 = rms_norm_src->nb[1] / ts0;
796    const int64_t s02 = rms_norm_src->nb[2] / ts0;
797    const int64_t s03 = rms_norm_src->nb[3] / ts0;
798
799    const size_t tdst = ggml_type_size(add_tensor->type);
800    GGML_ASSERT(add_tensor->nb[0] == tdst);
801    const int64_t d00 = add_tensor->nb[0] / tdst;
802    const int64_t d01 = add_tensor->nb[1] / tdst;
803    const int64_t d02 = add_tensor->nb[2] / tdst;
804    const int64_t d03 = add_tensor->nb[3] / tdst;
805
806    const size_t ts_mul = ggml_type_size(mul_src->type);
807    GGML_ASSERT(mul_src->nb[0] == ts_mul);
808    const int64_t mul_s01 = mul_src->nb[1] / ts_mul;
809    const int64_t mul_s02 = mul_src->nb[2] / ts_mul;
810    const int64_t mul_s03 = mul_src->nb[3] / ts_mul;
811    const int mul_nrows     = mul_src->ne[1];
812    const int mul_nchannels = mul_src->ne[2];
813    const int mul_nsamples  = mul_src->ne[3];
814
815    const size_t ts_add = ggml_type_size(add_src->type);
816    GGML_ASSERT(add_src->nb[0] == ts_add);
817    const int64_t add_s01 = add_src->nb[1] / ts_add;
818    const int64_t add_s02 = add_src->nb[2] / ts_add;
819    const int64_t add_s03 = add_src->nb[3] / ts_add;
820    const int add_nrows     = add_src->ne[1];
821    const int add_nchannels = add_src->ne[2];
822    const int add_nsamples  = add_src->ne[3];
823
824    rms_norm_mul_add_f32_sycl(src0_dd, mul_dd, add_dd, dst_dd, ne00, ne01, ne02, ne03,
825        s00, s01, s02, s03, d00, d01, d02, d03,
826        mul_s01, mul_s02, mul_s03, mul_nrows, mul_nchannels, mul_nsamples,
827        add_s01, add_s02, add_s03, add_nrows, add_nchannels, add_nsamples, eps, main_stream, ctx.device);
828}
829
830void ggml_sycl_op_rms_norm_back(ggml_backend_sycl_context & ctx, ggml_tensor * dst) {
831    scope_op_debug_print scope_dbg_print(__func__, dst, /*num_src=*/2);
832
833    GGML_ASSERT(dst->src[0]->type == GGML_TYPE_F32); // dz
834    GGML_ASSERT(dst->src[1]->type == GGML_TYPE_F32); // x
835    GGML_ASSERT(dst->type         == GGML_TYPE_F32);
836
837    float eps = 1e-5f;
838    std::memcpy(&eps, dst->op_params, sizeof(float));
839    if (!(eps > 0.0f) || !std::isfinite(eps)) eps = 1e-5f;
840
841    const float * g_base  = static_cast<const float *>(dst->src[0]->data); // dz
842    const float * x_base  = static_cast<const float *>(dst->src[1]->data); // x
843          float * dx_base = static_cast<      float *>(dst->data);
844
845    const int64_t D  = dst->ne[0];
846    const int64_t n1 = dst->ne[1], n2 = dst->ne[2], n3 = dst->ne[3]; (void) n3;
847    const int64_t N  = ggml_nrows(dst);
848    if (D == 0 || N == 0) return;
849
850    const ggml_tensor *G = dst->src[0];
851    const ggml_tensor *X = dst->src[1];
852    const int ts = (int) ggml_type_size(X->type);
853    GGML_ASSERT((size_t) X->nb[0]   == (size_t) ts);
854    GGML_ASSERT((size_t) G->nb[0]   == (size_t) ts);
855    GGML_ASSERT((size_t) dst->nb[0] == (size_t) ts);
856
857    const int64_t xs1 = X->nb[1] / ts, xs2 = X->nb[2] / ts, xs3 = X->nb[3] / ts;
858    const int64_t gs1 = G->nb[1] / ts, gs2 = G->nb[2] / ts, gs3 = G->nb[3] / ts;
859    const int64_t ds1 = dst->nb[1] / ts, ds2 = dst->nb[2] / ts, ds3 = dst->nb[3] / ts;
860
861    dpct::queue_ptr q = ctx.stream();
862
863    // work-group size: multiple of WARP_SIZE, capped by device and 256, and not larger than D
864    const int device_max_wg = ggml_sycl_info().max_work_group_sizes[ctx.device];
865    auto roundup = [](int v, int m) { return ((v + m - 1) / m) * m; };
866    int wg_cap = 256;
867    if (device_max_wg > 0) wg_cap = std::min(wg_cap, device_max_wg);
868    int WG = std::max(WARP_SIZE, std::min(roundup((int)std::min<int64_t>(D, wg_cap), WARP_SIZE), wg_cap));
869
870    // FP32 path: per-thread compensated accumulation + hierarchical reduction
871    q->submit([&](sycl::handler &cgh) {
872        const int nwarps_loc = std::max(1, WG / WARP_SIZE);
873        // store one partial value per warp (xx and xg) for cross-warp reduction
874        auto l_xx   = sycl::local_accessor<sycl::float2, 1>(sycl::range<1>(nwarps_loc), cgh);
875        auto l_xg   = sycl::local_accessor<sycl::float2, 1>(sycl::range<1>(nwarps_loc), cgh);
876
877        cgh.parallel_for(
878            sycl::nd_range<3>(sycl::range<3>(1, 1, N) * sycl::range<3>(1, 1, WG),
879                              sycl::range<3>(1, 1, WG)),
880            [=](sycl::nd_item<3> item_ct1) [[sycl::reqd_sub_group_size(WARP_SIZE)]] {
881                const int row = item_ct1.get_group(2);
882                const int tid = item_ct1.get_local_id(2);
883
884                const int64_t i1 = row % n1;
885                const int64_t i2 = (row / n1) % n2;
886                const int64_t i3 = row / (n1 * n2);
887
888                const float *__restrict x_row = x_base + i3 * xs3 + i2 * xs2 + i1 * xs1;
889                const float *__restrict g_row = g_base + i3 * gs3 + i2 * gs2 + i1 * gs1;
890                float *__restrict d_row       = dx_base + i3 * ds3 + i2 * ds2 + i1 * ds1;
891
892                // per-thread accumulation (compensated by default)
893                float sum_xx = 0.f, sum_xg = 0.f;
894#ifndef GGML_SYCL_RMS_BACK_FAST
895                float c_xx = 0.f, c_xg = 0.f;
896#endif
897                for (int64_t col = tid; col < D; col += WG) {
898                    const float xv = x_row[col];
899                    const float gv = g_row[col];
900#ifdef GGML_SYCL_RMS_BACK_FAST
901                    sum_xx += xv * xv;
902                    sum_xg += xv * gv;
903#else
904                    float y1 = xv * xv - c_xx;
905                    float t1 = sum_xx + y1;
906                    c_xx = (t1 - sum_xx) - y1;
907                    sum_xx = t1;
908
909                    float y2 = xv * gv - c_xg;
910                    float t2 = sum_xg + y2;
911                    c_xg = (t2 - sum_xg) - y2;
912                    sum_xg = t2;
913#endif
914                }
915
916                // warp-level reduction
917                sycl::float2 xx = sycl::float2(sum_xx,
918#ifndef GGML_SYCL_RMS_BACK_FAST
919                    c_xx
920#else
921                    0.f
922#endif
923                );
924                sycl::float2 xg = sycl::float2(sum_xg,
925#ifndef GGML_SYCL_RMS_BACK_FAST
926                    c_xg
927#else
928                    0.f
929#endif
930                );
931                xx = warp_reduce_sum(xx, item_ct1);
932                xg = warp_reduce_sum(xg, item_ct1);
933
934                // cross-warp reduction using local memory (single barrier)
935                const auto sub_group = item_ct1.get_sub_group();
936                const auto sg_id     = sub_group.get_group_linear_id();
937                const auto wi_in_sg  = sub_group.get_local_linear_id();
938                const int nthreads   = item_ct1.get_local_range(2);
939                const int nwarps     = nthreads / WARP_SIZE;
940
941                sycl::float2 xx_total = xx;
942                sycl::float2 xg_total = xg;
943                if (nwarps > 1) {
944                    if (wi_in_sg == 0) {
945                        l_xx[sg_id] = xx;
946                        l_xg[sg_id] = xg;
947                    }
948                    item_ct1.barrier(sycl::access::fence_space::local_space);
949
950                    if (sg_id == 0) {
951                        const unsigned wi_u = wi_in_sg;
952                        sycl::float2 xx_first = (wi_u < static_cast<unsigned>(nwarps)) ? l_xx[wi_u] : sycl::float2(0.f, 0.f);
953                        sycl::float2 xg_first = (wi_u < static_cast<unsigned>(nwarps)) ? l_xg[wi_u] : sycl::float2(0.f, 0.f);
954                        xx_total = warp_reduce_sum(xx_first, item_ct1);
955                        xg_total = warp_reduce_sum(xg_first, item_ct1);
956                    } else {
957                        // other subgroups keep their local totals; they'll be ignored
958                        xx_total = xx;
959                        xg_total = xg;
960                    }
961                    // ensure all threads see the first-subgroup result via broadcast below
962                }
963
964                // compute inv_r and coeff once per row and broadcast to the whole work-group
965                float inv_r = 0.f;
966                float coeff = 0.f;
967                if (tid == 0) {
968                    const float sum_xx_f  = xx_total.x() + xx_total.y();
969                    const float sum_xdz_f = xg_total.x() + xg_total.y();
970                    const float mean_eps  = sum_xx_f / (float) D + eps;
971                    const float sum_eps   = sum_xx_f + eps * (float) D;
972                    inv_r = sycl::rsqrt(mean_eps);
973                    coeff = -sum_xdz_f / sum_eps;
974                }
975                inv_r = sycl::group_broadcast(item_ct1.get_group(), inv_r);
976                coeff = sycl::group_broadcast(item_ct1.get_group(), coeff);
977
978                for (int64_t col = tid; col < D; col += WG) {
979                    d_row[col] = (g_row[col] + coeff * x_row[col]) * inv_r;
980                }
981            });
982    });
983
984}
985
986void ggml_sycl_op_l2_norm(ggml_backend_sycl_context& ctx, ggml_tensor* dst) {
987    const ggml_tensor * src0 = dst->src[0];
988    const float * src0_d = (const float *) src0->data;
989    float * dst_d = (float *) dst->data;
990    dpct::queue_ptr     stream = ctx.stream();
991
992    GGML_ASSERT(src0->type == GGML_TYPE_F32);
993    GGML_ASSERT( dst->type == GGML_TYPE_F32);
994
995    GGML_TENSOR_UNARY_OP_LOCALS;
996
997    float eps;
998    memcpy(&eps, dst->op_params, sizeof(float));
999    GGML_ASSERT(eps >= 0.0f);
1000
1001    const size_t ts0 = ggml_type_size(src0->type);
1002    const size_t tdst = ggml_type_size(dst->type);
1003    GGML_ASSERT(nb00 % ts0 == 0 && nb01 % ts0 == 0 && nb02 % ts0 == 0 && nb03 % ts0 == 0);
1004    GGML_ASSERT(nb0 % tdst == 0 && nb1 % tdst == 0 && nb2 % tdst == 0 && nb3 % tdst == 0);
1005    const int64_t ss0 = nb00 / ts0;
1006    const int64_t ss1 = nb01 / ts0;
1007    const int64_t ss2 = nb02 / ts0;
1008    const int64_t ss3 = nb03 / ts0;
1009    const int64_t ds0 = nb0 / tdst;
1010    const int64_t ds1 = nb1 / tdst;
1011    const int64_t ds2 = nb2 / tdst;
1012    const int64_t ds3 = nb3 / tdst;
1013
1014    /*support both WARP_SIZE or WARP_32_SIZE in code
1015      choose by hardware for better performance
1016    */
1017    l2_norm_f32_sycl<WARP_SIZE>(src0_d, dst_d, ne00, ne01, ne02, ne03,
1018            ss0, ss1, ss2, ss3, ds0, ds1, ds2, ds3, eps, stream, ctx.device);
1019}
1020
1021// nodes[0..count) are independent, same-shape, same-eps, same-nb L2_NORM ops validated
1022// by the caller; requires ncols < 1024 (the warp reduction path).
1023void ggml_sycl_l2_norm_batch(ggml_backend_sycl_context & ctx, ggml_tensor ** nodes, int count) {
1024    const ggml_tensor * s0 = nodes[0]->src[0];
1025    const int ncols       = (int) s0->ne[0];
1026    const int nrows_total = (int) ggml_nrows(s0);
1027    float eps;
1028    memcpy(&eps, nodes[0]->op_params, sizeof(float));
1029    GGML_ASSERT(eps >= 0.0f);
1030
1031    l2_batch_ptrs p{};
1032    for (int t = 0; t < count; ++t) {
1033        p.src[t] = (const float *) nodes[t]->src[0]->data;
1034        p.dst[t] = (float *) nodes[t]->data;
1035    }
1036
1037    const ggml_tensor * d0 = nodes[0];
1038    const size_t        ts = ggml_type_size(GGML_TYPE_F32);
1039    l2_batch_strides    st{};
1040    st.ne1 = (int) s0->ne[1];
1041    st.ne2 = (int) s0->ne[2];
1042    st.ss0 = s0->nb[0] / ts; st.ss1 = s0->nb[1] / ts; st.ss2 = s0->nb[2] / ts; st.ss3 = s0->nb[3] / ts;
1043    st.ds0 = d0->nb[0] / ts; st.ds1 = d0->nb[1] / ts; st.ds2 = d0->nb[2] / ts; st.ds3 = d0->nb[3] / ts;
1044
1045    l2_norm_f32_batch_sycl<WARP_SIZE>(p, st, count, ncols, nrows_total, eps, ctx.stream());
1046}