xref: /libCEED/backends/sycl-ref/ceed-sycl-vector.sycl.cpp (revision 9a25c3513918281076c2babe50808fbe5e3a546e)
1 // Copyright (c) 2017-2022, Lawrence Livermore National Security, LLC and other CEED contributors.
2 // All Rights Reserved. See the top-level LICENSE and NOTICE files for details.
3 //
4 // SPDX-License-Identifier: BSD-2-Clause
5 //
6 // This file is part of CEED:  http://github.com/ceed
7 
8 #include <ceed/backend.h>
9 #include <ceed/ceed.h>
10 
11 #include <cmath>
12 #include <string>
13 #include <sycl/sycl.hpp>
14 
15 #include "ceed-sycl-ref.hpp"
16 
17 //------------------------------------------------------------------------------
18 // Check if host/device sync is needed
19 //------------------------------------------------------------------------------
20 static inline int CeedVectorNeedSync_Sycl(const CeedVector vec, CeedMemType mem_type, bool *need_sync) {
21   bool             has_valid_array = false;
22   CeedVector_Sycl *impl;
23 
24   CeedCallBackend(CeedVectorGetData(vec, &impl));
25   CeedCallBackend(CeedVectorHasValidArray(vec, &has_valid_array));
26   switch (mem_type) {
27     case CEED_MEM_HOST:
28       *need_sync = has_valid_array && !impl->h_array;
29       break;
30     case CEED_MEM_DEVICE:
31       *need_sync = has_valid_array && !impl->d_array;
32       break;
33   }
34   return CEED_ERROR_SUCCESS;
35 }
36 
37 //------------------------------------------------------------------------------
38 // Sync host to device
39 //------------------------------------------------------------------------------
40 static inline int CeedVectorSyncH2D_Sycl(const CeedVector vec) {
41   Ceed             ceed;
42   Ceed_Sycl       *data;
43   CeedSize         length;
44   CeedVector_Sycl *impl;
45 
46   CeedCallBackend(CeedVectorGetCeed(vec, &ceed));
47   CeedCallBackend(CeedVectorGetData(vec, &impl));
48   CeedCallBackend(CeedGetData(ceed, &data));
49   CeedCheck(impl->h_array, ceed, CEED_ERROR_BACKEND, "No valid host data to sync to device");
50 
51   CeedCallBackend(CeedVectorGetLength(vec, &length));
52   if (impl->d_array_borrowed) {
53     impl->d_array = impl->d_array_borrowed;
54   } else if (impl->d_array_owned) {
55     impl->d_array = impl->d_array_owned;
56   } else {
57     CeedCallSycl(ceed, impl->d_array_owned = sycl::malloc_device<CeedScalar>(length, data->sycl_device, data->sycl_context));
58     impl->d_array = impl->d_array_owned;
59   }
60 
61   sycl::event e = data->sycl_queue.ext_oneapi_submit_barrier();
62   // Copy from host to device
63   sycl::event copy_event = data->sycl_queue.copy<CeedScalar>(impl->h_array, impl->d_array, length, {e});
64   // Wait for copy to finish and handle exceptions.
65   CeedCallSycl(ceed, copy_event.wait_and_throw());
66   return CEED_ERROR_SUCCESS;
67 }
68 
69 //------------------------------------------------------------------------------
70 // Sync device to host
71 //------------------------------------------------------------------------------
72 static inline int CeedVectorSyncD2H_Sycl(const CeedVector vec) {
73   Ceed             ceed;
74   Ceed_Sycl       *data;
75   CeedSize         length;
76   CeedVector_Sycl *impl;
77 
78   CeedCallBackend(CeedVectorGetCeed(vec, &ceed));
79   CeedCallBackend(CeedVectorGetData(vec, &impl));
80   CeedCallBackend(CeedGetData(ceed, &data));
81 
82   CeedCheck(impl->d_array, ceed, CEED_ERROR_BACKEND, "No valid device data to sync to host");
83 
84   CeedCallBackend(CeedVectorGetLength(vec, &length));
85   if (impl->h_array_borrowed) {
86     impl->h_array = impl->h_array_borrowed;
87   } else if (impl->h_array_owned) {
88     impl->h_array = impl->h_array_owned;
89   } else {
90     CeedCallBackend(CeedCalloc(length, &impl->h_array_owned));
91     impl->h_array = impl->h_array_owned;
92   }
93 
94   // Order queue
95   sycl::event e = data->sycl_queue.ext_oneapi_submit_barrier();
96   // Copy from device to host
97   sycl::event copy_event = data->sycl_queue.copy<CeedScalar>(impl->d_array, impl->h_array, length, {e});
98   // Wait for copy to finish and handle exceptions.
99   CeedCallSycl(ceed, copy_event.wait_and_throw());
100   return CEED_ERROR_SUCCESS;
101 }
102 
103 //------------------------------------------------------------------------------
104 // Sync arrays
105 //------------------------------------------------------------------------------
106 static int CeedVectorSyncArray_Sycl(const CeedVector vec, CeedMemType mem_type) {
107   bool need_sync = false;
108 
109   // Check whether device/host sync is needed
110   CeedCallBackend(CeedVectorNeedSync_Sycl(vec, mem_type, &need_sync));
111   if (!need_sync) return CEED_ERROR_SUCCESS;
112 
113   switch (mem_type) {
114     case CEED_MEM_HOST:
115       return CeedVectorSyncD2H_Sycl(vec);
116     case CEED_MEM_DEVICE:
117       return CeedVectorSyncH2D_Sycl(vec);
118   }
119   return CEED_ERROR_UNSUPPORTED;
120 }
121 
122 //------------------------------------------------------------------------------
123 // Set all pointers as invalid
124 //------------------------------------------------------------------------------
125 static inline int CeedVectorSetAllInvalid_Sycl(const CeedVector vec) {
126   CeedVector_Sycl *impl;
127 
128   CeedCallBackend(CeedVectorGetData(vec, &impl));
129   impl->h_array = NULL;
130   impl->d_array = NULL;
131   return CEED_ERROR_SUCCESS;
132 }
133 
134 //------------------------------------------------------------------------------
135 // Check if CeedVector has any valid pointer
136 //------------------------------------------------------------------------------
137 static inline int CeedVectorHasValidArray_Sycl(const CeedVector vec, bool *has_valid_array) {
138   CeedVector_Sycl *impl;
139 
140   CeedCallBackend(CeedVectorGetData(vec, &impl));
141   *has_valid_array = impl->h_array || impl->d_array;
142   return CEED_ERROR_SUCCESS;
143 }
144 
145 //------------------------------------------------------------------------------
146 // Check if has array of given type
147 //------------------------------------------------------------------------------
148 static inline int CeedVectorHasArrayOfType_Sycl(const CeedVector vec, CeedMemType mem_type, bool *has_array_of_type) {
149   CeedVector_Sycl *impl;
150 
151   CeedCallBackend(CeedVectorGetData(vec, &impl));
152   switch (mem_type) {
153     case CEED_MEM_HOST:
154       *has_array_of_type = impl->h_array_borrowed || impl->h_array_owned;
155       break;
156     case CEED_MEM_DEVICE:
157       *has_array_of_type = impl->d_array_borrowed || impl->d_array_owned;
158       break;
159   }
160   return CEED_ERROR_SUCCESS;
161 }
162 
163 //------------------------------------------------------------------------------
164 // Check if has borrowed array of given type
165 //------------------------------------------------------------------------------
166 static inline int CeedVectorHasBorrowedArrayOfType_Sycl(const CeedVector vec, CeedMemType mem_type, bool *has_borrowed_array_of_type) {
167   CeedVector_Sycl *impl;
168 
169   CeedCallBackend(CeedVectorGetData(vec, &impl));
170   switch (mem_type) {
171     case CEED_MEM_HOST:
172       *has_borrowed_array_of_type = impl->h_array_borrowed;
173       break;
174     case CEED_MEM_DEVICE:
175       *has_borrowed_array_of_type = impl->d_array_borrowed;
176       break;
177   }
178   return CEED_ERROR_SUCCESS;
179 }
180 
181 //------------------------------------------------------------------------------
182 // Set array from host
183 //------------------------------------------------------------------------------
184 static int CeedVectorSetArrayHost_Sycl(const CeedVector vec, const CeedCopyMode copy_mode, CeedScalar *array) {
185   CeedVector_Sycl *impl;
186 
187   CeedCallBackend(CeedVectorGetData(vec, &impl));
188   switch (copy_mode) {
189     case CEED_COPY_VALUES: {
190       if (!impl->h_array_owned) {
191         CeedSize length;
192 
193         CeedCallBackend(CeedVectorGetLength(vec, &length));
194         CeedCallBackend(CeedMalloc(length, &impl->h_array_owned));
195       }
196       impl->h_array_borrowed = NULL;
197       impl->h_array          = impl->h_array_owned;
198       if (array) {
199         CeedSize length;
200         size_t   bytes;
201 
202         CeedCallBackend(CeedVectorGetLength(vec, &length));
203         bytes = length * sizeof(CeedScalar);
204         memcpy(impl->h_array, array, bytes);
205       }
206     } break;
207     case CEED_OWN_POINTER:
208       CeedCallBackend(CeedFree(&impl->h_array_owned));
209       impl->h_array_owned    = array;
210       impl->h_array_borrowed = NULL;
211       impl->h_array          = array;
212       break;
213     case CEED_USE_POINTER:
214       CeedCallBackend(CeedFree(&impl->h_array_owned));
215       impl->h_array_borrowed = array;
216       impl->h_array          = array;
217       break;
218   }
219   return CEED_ERROR_SUCCESS;
220 }
221 
222 //------------------------------------------------------------------------------
223 // Set array from device
224 //------------------------------------------------------------------------------
225 static int CeedVectorSetArrayDevice_Sycl(const CeedVector vec, const CeedCopyMode copy_mode, CeedScalar *array) {
226   Ceed             ceed;
227   Ceed_Sycl       *data;
228   CeedVector_Sycl *impl;
229 
230   CeedCallBackend(CeedVectorGetCeed(vec, &ceed));
231   CeedCallBackend(CeedVectorGetData(vec, &impl));
232   CeedCallBackend(CeedGetData(ceed, &data));
233 
234   // Order queue
235   sycl::event e = data->sycl_queue.ext_oneapi_submit_barrier();
236 
237   switch (copy_mode) {
238     case CEED_COPY_VALUES: {
239       CeedSize length;
240 
241       CeedCallBackend(CeedVectorGetLength(vec, &length));
242       if (!impl->d_array_owned) {
243         CeedCallSycl(ceed, impl->d_array_owned = sycl::malloc_device<CeedScalar>(length, data->sycl_device, data->sycl_context));
244       }
245       impl->d_array = impl->d_array_owned;
246       if (array) {
247         sycl::event copy_event = data->sycl_queue.copy<CeedScalar>(array, impl->d_array, length, {e});
248         // Wait for copy to finish and handle exceptions.
249         CeedCallSycl(ceed, copy_event.wait_and_throw());
250       }
251     } break;
252     case CEED_OWN_POINTER:
253       if (impl->d_array_owned) {
254         // Wait for all work to finish before freeing memory
255         CeedCallSycl(ceed, data->sycl_queue.wait_and_throw());
256         CeedCallSycl(ceed, sycl::free(impl->d_array_owned, data->sycl_context));
257       }
258       impl->d_array_owned    = array;
259       impl->d_array_borrowed = NULL;
260       impl->d_array          = array;
261       break;
262     case CEED_USE_POINTER:
263       if (impl->d_array_owned) {
264         // Wait for all work to finish before freeing memory
265         CeedCallSycl(ceed, data->sycl_queue.wait_and_throw());
266         CeedCallSycl(ceed, sycl::free(impl->d_array_owned, data->sycl_context));
267       }
268       impl->d_array_owned    = NULL;
269       impl->d_array_borrowed = array;
270       impl->d_array          = array;
271       break;
272   }
273   return CEED_ERROR_SUCCESS;
274 }
275 
276 //------------------------------------------------------------------------------
277 // Set the array used by a vector,
278 //   freeing any previously allocated array if applicable
279 //------------------------------------------------------------------------------
280 static int CeedVectorSetArray_Sycl(const CeedVector vec, const CeedMemType mem_type, const CeedCopyMode copy_mode, CeedScalar *array) {
281   Ceed             ceed;
282   CeedVector_Sycl *impl;
283 
284   CeedCallBackend(CeedVectorGetCeed(vec, &ceed));
285   CeedCallBackend(CeedVectorGetData(vec, &impl));
286 
287   CeedCallBackend(CeedVectorSetAllInvalid_Sycl(vec));
288   switch (mem_type) {
289     case CEED_MEM_HOST:
290       return CeedVectorSetArrayHost_Sycl(vec, copy_mode, array);
291     case CEED_MEM_DEVICE:
292       return CeedVectorSetArrayDevice_Sycl(vec, copy_mode, array);
293   }
294   return CEED_ERROR_UNSUPPORTED;
295 }
296 
297 //------------------------------------------------------------------------------
298 // Set host array to value
299 //------------------------------------------------------------------------------
300 static int CeedHostSetValue_Sycl(CeedScalar *h_array, CeedSize length, CeedScalar val) {
301   for (CeedSize i = 0; i < length; i++) h_array[i] = val;
302   return CEED_ERROR_SUCCESS;
303 }
304 
305 //------------------------------------------------------------------------------
306 // Set device array to value
307 //------------------------------------------------------------------------------
308 static int CeedDeviceSetValue_Sycl(sycl::queue &sycl_queue, CeedScalar *d_array, CeedSize length, CeedScalar val) {
309   // Order queue
310   sycl::event e = sycl_queue.ext_oneapi_submit_barrier();
311   sycl_queue.fill(d_array, val, length, {e});
312   return CEED_ERROR_SUCCESS;
313 }
314 
315 //------------------------------------------------------------------------------
316 // Set a vector to a value,
317 //------------------------------------------------------------------------------
318 static int CeedVectorSetValue_Sycl(CeedVector vec, CeedScalar val) {
319   Ceed             ceed;
320   Ceed_Sycl       *data;
321   CeedSize         length;
322   CeedVector_Sycl *impl;
323 
324   CeedCallBackend(CeedVectorGetCeed(vec, &ceed));
325   CeedCallBackend(CeedVectorGetData(vec, &impl));
326   CeedCallBackend(CeedVectorGetLength(vec, &length));
327   CeedCallBackend(CeedGetData(ceed, &data));
328 
329   // Set value for synced device/host array
330   if (!impl->d_array && !impl->h_array) {
331     if (impl->d_array_borrowed) {
332       impl->d_array = impl->d_array_borrowed;
333     } else if (impl->h_array_borrowed) {
334       impl->h_array = impl->h_array_borrowed;
335     } else if (impl->d_array_owned) {
336       impl->d_array = impl->d_array_owned;
337     } else if (impl->h_array_owned) {
338       impl->h_array = impl->h_array_owned;
339     } else {
340       CeedCallBackend(CeedVectorSetArray(vec, CEED_MEM_DEVICE, CEED_COPY_VALUES, NULL));
341     }
342   }
343   if (impl->d_array) {
344     CeedCallBackend(CeedDeviceSetValue_Sycl(data->sycl_queue, impl->d_array, length, val));
345     impl->h_array = NULL;
346   }
347   if (impl->h_array) {
348     CeedCallBackend(CeedHostSetValue_Sycl(impl->h_array, length, val));
349     impl->d_array = NULL;
350   }
351   return CEED_ERROR_SUCCESS;
352 }
353 
354 //------------------------------------------------------------------------------
355 // Vector Take Array
356 //------------------------------------------------------------------------------
357 static int CeedVectorTakeArray_Sycl(CeedVector vec, CeedMemType mem_type, CeedScalar **array) {
358   Ceed             ceed;
359   Ceed_Sycl       *data;
360   CeedVector_Sycl *impl;
361 
362   CeedCallBackend(CeedVectorGetCeed(vec, &ceed));
363   CeedCallBackend(CeedVectorGetData(vec, &impl));
364   CeedCallBackend(CeedGetData(ceed, &data));
365 
366   // Order queue
367   data->sycl_queue.ext_oneapi_submit_barrier();
368 
369   // Sync array to requested mem_type
370   CeedCallBackend(CeedVectorSyncArray(vec, mem_type));
371 
372   // Update pointer
373   switch (mem_type) {
374     case CEED_MEM_HOST:
375       (*array)               = impl->h_array_borrowed;
376       impl->h_array_borrowed = NULL;
377       impl->h_array          = NULL;
378       break;
379     case CEED_MEM_DEVICE:
380       (*array)               = impl->d_array_borrowed;
381       impl->d_array_borrowed = NULL;
382       impl->d_array          = NULL;
383       break;
384   }
385   return CEED_ERROR_SUCCESS;
386 }
387 
388 //------------------------------------------------------------------------------
389 // Core logic for array syncronization for GetArray.
390 //   If a different memory type is most up to date, this will perform a copy
391 //------------------------------------------------------------------------------
392 static int CeedVectorGetArrayCore_Sycl(const CeedVector vec, const CeedMemType mem_type, CeedScalar **array) {
393   Ceed             ceed;
394   CeedVector_Sycl *impl;
395 
396   CeedCallBackend(CeedVectorGetCeed(vec, &ceed));
397   CeedCallBackend(CeedVectorGetData(vec, &impl));
398 
399   // Sync array to requested mem_type
400   CeedCallBackend(CeedVectorSyncArray(vec, mem_type));
401 
402   // Update pointer
403   switch (mem_type) {
404     case CEED_MEM_HOST:
405       *array = impl->h_array;
406       break;
407     case CEED_MEM_DEVICE:
408       *array = impl->d_array;
409       break;
410   }
411   return CEED_ERROR_SUCCESS;
412 }
413 
414 //------------------------------------------------------------------------------
415 // Get read-only access to a vector via the specified mem_type
416 //------------------------------------------------------------------------------
417 static int CeedVectorGetArrayRead_Sycl(const CeedVector vec, const CeedMemType mem_type, const CeedScalar **array) {
418   return CeedVectorGetArrayCore_Sycl(vec, mem_type, (CeedScalar **)array);
419 }
420 
421 //------------------------------------------------------------------------------
422 // Get read/write access to a vector via the specified mem_type
423 //------------------------------------------------------------------------------
424 static int CeedVectorGetArray_Sycl(const CeedVector vec, const CeedMemType mem_type, CeedScalar **array) {
425   CeedVector_Sycl *impl;
426 
427   CeedCallBackend(CeedVectorGetData(vec, &impl));
428   CeedCallBackend(CeedVectorGetArrayCore_Sycl(vec, mem_type, array));
429   CeedCallBackend(CeedVectorSetAllInvalid_Sycl(vec));
430   switch (mem_type) {
431     case CEED_MEM_HOST:
432       impl->h_array = *array;
433       break;
434     case CEED_MEM_DEVICE:
435       impl->d_array = *array;
436       break;
437   }
438   return CEED_ERROR_SUCCESS;
439 }
440 
441 //------------------------------------------------------------------------------
442 // Get write access to a vector via the specified mem_type
443 //------------------------------------------------------------------------------
444 static int CeedVectorGetArrayWrite_Sycl(const CeedVector vec, const CeedMemType mem_type, CeedScalar **array) {
445   bool             has_array_of_type = true;
446   CeedVector_Sycl *impl;
447 
448   CeedCallBackend(CeedVectorGetData(vec, &impl));
449   CeedCallBackend(CeedVectorHasArrayOfType_Sycl(vec, mem_type, &has_array_of_type));
450   if (!has_array_of_type) {
451     // Allocate if array is not yet allocated
452     CeedCallBackend(CeedVectorSetArray(vec, mem_type, CEED_COPY_VALUES, NULL));
453   } else {
454     // Select dirty array
455     switch (mem_type) {
456       case CEED_MEM_HOST:
457         if (impl->h_array_borrowed) impl->h_array = impl->h_array_borrowed;
458         else impl->h_array = impl->h_array_owned;
459         break;
460       case CEED_MEM_DEVICE:
461         if (impl->d_array_borrowed) impl->d_array = impl->d_array_borrowed;
462         else impl->d_array = impl->d_array_owned;
463     }
464   }
465   return CeedVectorGetArray_Sycl(vec, mem_type, array);
466 }
467 
468 //------------------------------------------------------------------------------
469 // Get the norm of a CeedVector
470 //------------------------------------------------------------------------------
471 static int CeedVectorNorm_Sycl(CeedVector vec, CeedNormType type, CeedScalar *norm) {
472   Ceed              ceed;
473   Ceed_Sycl        *data;
474   CeedSize          length;
475   const CeedScalar *d_array;
476   CeedVector_Sycl  *impl;
477 
478   CeedCallBackend(CeedVectorGetCeed(vec, &ceed));
479   CeedCallBackend(CeedVectorGetData(vec, &impl));
480   CeedCallBackend(CeedVectorGetLength(vec, &length));
481   CeedCallBackend(CeedGetData(ceed, &data));
482 
483   // Compute norm
484   CeedCallBackend(CeedVectorGetArrayRead(vec, CEED_MEM_DEVICE, &d_array));
485   switch (type) {
486     case CEED_NORM_1: {
487       // Order queue
488       sycl::event e            = data->sycl_queue.ext_oneapi_submit_barrier();
489       auto        sumReduction = sycl::reduction(impl->reduction_norm, sycl::plus<>(), {sycl::property::reduction::initialize_to_identity{}});
490       data->sycl_queue.parallel_for(length, {e}, sumReduction, [=](sycl::id<1> i, auto &sum) { sum += abs(d_array[i]); }).wait_and_throw();
491     } break;
492     case CEED_NORM_2: {
493       // Order queue
494       sycl::event e            = data->sycl_queue.ext_oneapi_submit_barrier();
495       auto        sumReduction = sycl::reduction(impl->reduction_norm, sycl::plus<>(), {sycl::property::reduction::initialize_to_identity{}});
496       data->sycl_queue.parallel_for(length, {e}, sumReduction, [=](sycl::id<1> i, auto &sum) { sum += (d_array[i] * d_array[i]); }).wait_and_throw();
497     } break;
498     case CEED_NORM_MAX: {
499       // Order queue
500       sycl::event e            = data->sycl_queue.ext_oneapi_submit_barrier();
501       auto        maxReduction = sycl::reduction(impl->reduction_norm, sycl::maximum<>(), {sycl::property::reduction::initialize_to_identity{}});
502       data->sycl_queue.parallel_for(length, {e}, maxReduction, [=](sycl::id<1> i, auto &max) { max.combine(abs(d_array[i])); }).wait_and_throw();
503     } break;
504   }
505   // L2 norm - square root over reduced value
506   if (type == CEED_NORM_2) *norm = sqrt(*impl->reduction_norm);
507   else *norm = *impl->reduction_norm;
508   CeedCallBackend(CeedVectorRestoreArrayRead(vec, &d_array));
509   return CEED_ERROR_SUCCESS;
510 }
511 
512 //------------------------------------------------------------------------------
513 // Take reciprocal of a vector on host
514 //------------------------------------------------------------------------------
515 static int CeedHostReciprocal_Sycl(CeedScalar *h_array, CeedSize length) {
516   for (CeedSize i = 0; i < length; i++) {
517     if (std::fabs(h_array[i]) > CEED_EPSILON) h_array[i] = 1. / h_array[i];
518   }
519   return CEED_ERROR_SUCCESS;
520 }
521 
522 //------------------------------------------------------------------------------
523 // Take reciprocal of a vector on device
524 //------------------------------------------------------------------------------
525 static int CeedDeviceReciprocal_Sycl(sycl::queue &sycl_queue, CeedScalar *d_array, CeedSize length) {
526   // Order queue
527   sycl::event e = sycl_queue.ext_oneapi_submit_barrier();
528   sycl_queue.parallel_for(length, {e}, [=](sycl::id<1> i) {
529     if (std::fabs(d_array[i]) > CEED_EPSILON) d_array[i] = 1. / d_array[i];
530   });
531   return CEED_ERROR_SUCCESS;
532 }
533 
534 //------------------------------------------------------------------------------
535 // Take reciprocal of a vector
536 //------------------------------------------------------------------------------
537 static int CeedVectorReciprocal_Sycl(CeedVector vec) {
538   Ceed             ceed;
539   Ceed_Sycl       *data;
540   CeedSize         length;
541   CeedVector_Sycl *impl;
542 
543   CeedCallBackend(CeedVectorGetCeed(vec, &ceed));
544   CeedCallBackend(CeedVectorGetData(vec, &impl));
545   CeedCallBackend(CeedVectorGetLength(vec, &length));
546   CeedCallBackend(CeedGetData(ceed, &data));
547 
548   // Set value for synced device/host array
549   if (impl->d_array) CeedCallBackend(CeedDeviceReciprocal_Sycl(data->sycl_queue, impl->d_array, length));
550   if (impl->h_array) CeedCallBackend(CeedHostReciprocal_Sycl(impl->h_array, length));
551   return CEED_ERROR_SUCCESS;
552 }
553 
554 //------------------------------------------------------------------------------
555 // Compute x = alpha x on the host
556 //------------------------------------------------------------------------------
557 static int CeedHostScale_Sycl(CeedScalar *x_array, CeedScalar alpha, CeedSize length) {
558   for (CeedSize i = 0; i < length; i++) x_array[i] *= alpha;
559   return CEED_ERROR_SUCCESS;
560 }
561 
562 //------------------------------------------------------------------------------
563 // Compute x = alpha x on device
564 //------------------------------------------------------------------------------
565 static int CeedDeviceScale_Sycl(sycl::queue &sycl_queue, CeedScalar *x_array, CeedScalar alpha, CeedSize length) {
566   // Order queue
567   sycl::event e = sycl_queue.ext_oneapi_submit_barrier();
568   sycl_queue.parallel_for(length, {e}, [=](sycl::id<1> i) { x_array[i] *= alpha; });
569   return CEED_ERROR_SUCCESS;
570 }
571 
572 //------------------------------------------------------------------------------
573 // Compute x = alpha x
574 //------------------------------------------------------------------------------
575 static int CeedVectorScale_Sycl(CeedVector x, CeedScalar alpha) {
576   Ceed             ceed;
577   Ceed_Sycl       *data;
578   CeedSize         length;
579   CeedVector_Sycl *x_impl;
580 
581   CeedCallBackend(CeedVectorGetCeed(x, &ceed));
582   CeedCallBackend(CeedVectorGetData(x, &x_impl));
583   CeedCallBackend(CeedVectorGetLength(x, &length));
584   CeedCallBackend(CeedGetData(ceed, &data));
585 
586   // Set value for synced device/host array
587   if (x_impl->d_array) CeedCallBackend(CeedDeviceScale_Sycl(data->sycl_queue, x_impl->d_array, alpha, length));
588   if (x_impl->h_array) CeedCallBackend(CeedHostScale_Sycl(x_impl->h_array, alpha, length));
589   return CEED_ERROR_SUCCESS;
590 }
591 
592 //------------------------------------------------------------------------------
593 // Compute y = alpha x + y on the host
594 //------------------------------------------------------------------------------
595 static int CeedHostAXPY_Sycl(CeedScalar *y_array, CeedScalar alpha, CeedScalar *x_array, CeedSize length) {
596   for (CeedSize i = 0; i < length; i++) y_array[i] += alpha * x_array[i];
597   return CEED_ERROR_SUCCESS;
598 }
599 
600 //------------------------------------------------------------------------------
601 // Compute y = alpha x + y on device
602 //------------------------------------------------------------------------------
603 static int CeedDeviceAXPY_Sycl(sycl::queue &sycl_queue, CeedScalar *y_array, CeedScalar alpha, CeedScalar *x_array, CeedSize length) {
604   // Order queue
605   sycl::event e = sycl_queue.ext_oneapi_submit_barrier();
606   sycl_queue.parallel_for(length, {e}, [=](sycl::id<1> i) { y_array[i] += alpha * x_array[i]; });
607   return CEED_ERROR_SUCCESS;
608 }
609 
610 //------------------------------------------------------------------------------
611 // Compute y = alpha x + y
612 //------------------------------------------------------------------------------
613 static int CeedVectorAXPY_Sycl(CeedVector y, CeedScalar alpha, CeedVector x) {
614   Ceed             ceed;
615   Ceed_Sycl       *data;
616   CeedSize         length;
617   CeedVector_Sycl *y_impl, *x_impl;
618 
619   CeedCallBackend(CeedVectorGetCeed(y, &ceed));
620   CeedCallBackend(CeedVectorGetData(y, &y_impl));
621   CeedCallBackend(CeedVectorGetData(x, &x_impl));
622   CeedCallBackend(CeedVectorGetLength(y, &length));
623   CeedCallBackend(CeedGetData(ceed, &data));
624 
625   // Set value for synced device/host array
626   if (y_impl->d_array) {
627     CeedCallBackend(CeedVectorSyncArray(x, CEED_MEM_DEVICE));
628     CeedCallBackend(CeedDeviceAXPY_Sycl(data->sycl_queue, y_impl->d_array, alpha, x_impl->d_array, length));
629   }
630   if (y_impl->h_array) {
631     CeedCallBackend(CeedVectorSyncArray(x, CEED_MEM_HOST));
632     CeedCallBackend(CeedHostAXPY_Sycl(y_impl->h_array, alpha, x_impl->h_array, length));
633   }
634   return CEED_ERROR_SUCCESS;
635 }
636 
637 //------------------------------------------------------------------------------
638 // Compute the pointwise multiplication w = x .* y on the host
639 //------------------------------------------------------------------------------
640 static int CeedHostPointwiseMult_Sycl(CeedScalar *w_array, CeedScalar *x_array, CeedScalar *y_array, CeedSize length) {
641   for (CeedSize i = 0; i < length; i++) w_array[i] = x_array[i] * y_array[i];
642   return CEED_ERROR_SUCCESS;
643 }
644 
645 //------------------------------------------------------------------------------
646 // Compute the pointwise multiplication w = x .* y on device (impl in .cu file)
647 //------------------------------------------------------------------------------
648 static int CeedDevicePointwiseMult_Sycl(sycl::queue &sycl_queue, CeedScalar *w_array, CeedScalar *x_array, CeedScalar *y_array, CeedSize length) {
649   // Order queue
650   sycl::event e = sycl_queue.ext_oneapi_submit_barrier();
651   sycl_queue.parallel_for(length, {e}, [=](sycl::id<1> i) { w_array[i] = x_array[i] * y_array[i]; });
652   return CEED_ERROR_SUCCESS;
653 }
654 
655 //------------------------------------------------------------------------------
656 // Compute the pointwise multiplication w = x .* y
657 //------------------------------------------------------------------------------
658 static int CeedVectorPointwiseMult_Sycl(CeedVector w, CeedVector x, CeedVector y) {
659   Ceed             ceed;
660   Ceed_Sycl       *data;
661   CeedSize         length;
662   CeedVector_Sycl *w_impl, *x_impl, *y_impl;
663 
664   CeedCallBackend(CeedVectorGetCeed(w, &ceed));
665   CeedCallBackend(CeedVectorGetData(w, &w_impl));
666   CeedCallBackend(CeedVectorGetData(x, &x_impl));
667   CeedCallBackend(CeedVectorGetData(y, &y_impl));
668   CeedCallBackend(CeedVectorGetLength(w, &length));
669   CeedCallBackend(CeedGetData(ceed, &data));
670 
671   // Set value for synced device/host array
672   if (!w_impl->d_array && !w_impl->h_array) {
673     CeedCallBackend(CeedVectorSetValue(w, 0.0));
674   }
675   if (w_impl->d_array) {
676     CeedCallBackend(CeedVectorSyncArray(x, CEED_MEM_DEVICE));
677     CeedCallBackend(CeedVectorSyncArray(y, CEED_MEM_DEVICE));
678     CeedCallBackend(CeedDevicePointwiseMult_Sycl(data->sycl_queue, w_impl->d_array, x_impl->d_array, y_impl->d_array, length));
679   }
680   if (w_impl->h_array) {
681     CeedCallBackend(CeedVectorSyncArray(x, CEED_MEM_HOST));
682     CeedCallBackend(CeedVectorSyncArray(y, CEED_MEM_HOST));
683     CeedCallBackend(CeedHostPointwiseMult_Sycl(w_impl->h_array, x_impl->h_array, y_impl->h_array, length));
684   }
685   return CEED_ERROR_SUCCESS;
686 }
687 
688 //------------------------------------------------------------------------------
689 // Destroy the vector
690 //------------------------------------------------------------------------------
691 static int CeedVectorDestroy_Sycl(const CeedVector vec) {
692   Ceed             ceed;
693   Ceed_Sycl       *data;
694   CeedVector_Sycl *impl;
695 
696   CeedCallBackend(CeedVectorGetCeed(vec, &ceed));
697   CeedCallBackend(CeedVectorGetData(vec, &impl));
698   CeedCallBackend(CeedGetData(ceed, &data));
699 
700   // Wait for all work to finish before freeing memory
701   CeedCallSycl(ceed, data->sycl_queue.wait_and_throw());
702   CeedCallSycl(ceed, sycl::free(impl->d_array_owned, data->sycl_context));
703   CeedCallSycl(ceed, sycl::free(impl->reduction_norm, data->sycl_context));
704 
705   CeedCallBackend(CeedFree(&impl->h_array_owned));
706   CeedCallBackend(CeedFree(&impl));
707   return CEED_ERROR_SUCCESS;
708 }
709 
710 //------------------------------------------------------------------------------
711 // Create a vector of the specified length (does not allocate memory)
712 //------------------------------------------------------------------------------
713 int CeedVectorCreate_Sycl(CeedSize n, CeedVector vec) {
714   Ceed             ceed;
715   Ceed_Sycl       *data;
716   CeedVector_Sycl *impl;
717 
718   CeedCallBackend(CeedVectorGetCeed(vec, &ceed));
719   CeedCallBackend(CeedGetData(ceed, &data));
720   CeedCallBackend(CeedCalloc(1, &impl));
721   CeedCallSycl(ceed, impl->reduction_norm = sycl::malloc_host<CeedScalar>(1, data->sycl_context));
722   CeedCallBackend(CeedSetBackendFunctionCpp(ceed, "Vector", vec, "HasValidArray", CeedVectorHasValidArray_Sycl));
723   CeedCallBackend(CeedSetBackendFunctionCpp(ceed, "Vector", vec, "HasBorrowedArrayOfType", CeedVectorHasBorrowedArrayOfType_Sycl));
724   CeedCallBackend(CeedSetBackendFunctionCpp(ceed, "Vector", vec, "SetArray", CeedVectorSetArray_Sycl));
725   CeedCallBackend(CeedSetBackendFunctionCpp(ceed, "Vector", vec, "TakeArray", CeedVectorTakeArray_Sycl));
726   CeedCallBackend(CeedSetBackendFunctionCpp(ceed, "Vector", vec, "SetValue", CeedVectorSetValue_Sycl));
727   CeedCallBackend(CeedSetBackendFunctionCpp(ceed, "Vector", vec, "SyncArray", CeedVectorSyncArray_Sycl));
728   CeedCallBackend(CeedSetBackendFunctionCpp(ceed, "Vector", vec, "GetArray", CeedVectorGetArray_Sycl));
729   CeedCallBackend(CeedSetBackendFunctionCpp(ceed, "Vector", vec, "GetArrayRead", CeedVectorGetArrayRead_Sycl));
730   CeedCallBackend(CeedSetBackendFunctionCpp(ceed, "Vector", vec, "GetArrayWrite", CeedVectorGetArrayWrite_Sycl));
731   CeedCallBackend(CeedSetBackendFunctionCpp(ceed, "Vector", vec, "Norm", CeedVectorNorm_Sycl));
732   CeedCallBackend(CeedSetBackendFunctionCpp(ceed, "Vector", vec, "Reciprocal", CeedVectorReciprocal_Sycl));
733   CeedCallBackend(CeedSetBackendFunctionCpp(ceed, "Vector", vec, "AXPY", CeedVectorAXPY_Sycl));
734   CeedCallBackend(CeedSetBackendFunctionCpp(ceed, "Vector", vec, "Scale", CeedVectorScale_Sycl));
735   CeedCallBackend(CeedSetBackendFunctionCpp(ceed, "Vector", vec, "PointwiseMult", CeedVectorPointwiseMult_Sycl));
736   CeedCallBackend(CeedSetBackendFunctionCpp(ceed, "Vector", vec, "Destroy", CeedVectorDestroy_Sycl));
737   CeedCallBackend(CeedVectorSetData(vec, impl));
738   return CEED_ERROR_SUCCESS;
739 }
740