Alien  1.3.0
User documentation
Loading...
Searching...
No Matches
SYCLInternalLinearAlgebra.cc
1// -*- tab-width: 2; indent-tabs-mode: nil; coding: utf-8-with-signature -*-
2//-----------------------------------------------------------------------------
3// Copyright 2000-2026 CEA (www.cea.fr) IFPEN (www.ifpenergiesnouvelles.com)
4// See the top-level COPYRIGHT file for details.
5// SPDX-License-Identifier: Apache-2.0
6//-----------------------------------------------------------------------------
7
8
9#include "arccore/message_passing/ITypeDispatcher.h"
10#include "arccore/message_passing/Request.h"
11#include "arccore/message_passing/IStat.h"
12
13#include <alien/utils/Precomp.h>
14
15#include <arccore/base/NotImplementedException.h>
16#include <arccore/base/TraceInfo.h>
17
20
21#include "alien/kernels/sycl/data/SYCLEnv.h"
22#include "alien/kernels/sycl/data/SYCLEnvInternal.h"
23
24#include "alien/kernels/sycl/data/SYCLVectorInternal.h"
25#include "alien/kernels/sycl/data/SYCLVector.h"
26
27#include <alien/kernels/sycl/data/SYCLBEllPackInternal.h>
28#include <alien/kernels/sycl/data/SYCLBEllPackMatrix.h>
29
30#include "alien/kernels/sycl/algebra/SYCLKernelInternal.h"
31#include "alien/kernels/sycl/algebra/SYCLInternalLinearAlgebra.h"
32
33#include "alien/kernels/sycl/algebra/SYCLBEllPackMatrixMult.h"
34
35/*---------------------------------------------------------------------------*/
36
37namespace Alien
38{
39
40using namespace Arccore;
41
42template class ALIEN_EXPORT LinearAlgebra<BackEnd::tag::sycl>;
43template class ALIEN_EXPORT LinearAlgebraExpr<BackEnd::tag::sycl>;
44
46SYCLInternalLinearAlgebraFactory()
47{
48 return new SYCLInternalLinearAlgebra();
49}
50
52SYCLInternalLinearAlgebraExprFactory()
53{
55}
56
57/*---------------------------------------------------------------------------*/
58
59namespace Internal = SYCLInternal;
60
61/*---------------------------------------------------------------------------*/
62
63SYCLInternalLinearAlgebra::SYCLInternalLinearAlgebra()
65{
66 m_internal.reset(new SYCLInternal::KernelInternal());
67}
68
69/*---------------------------------------------------------------------------*/
70
71SYCLInternalLinearAlgebra::~SYCLInternalLinearAlgebra()
72{
73#ifdef ALIEN_USE_PERF_TIMER
74 m_timer.printInfo("SYCL-ALGEBRA");
75#endif
76}
77
78void SYCLInternalLinearAlgebra::setDotAlgo(int dot_algo)
79{
80 m_internal->setDotAlgo(dot_algo);
81}
82
83/*---------------------------------------------------------------------------*/
84/*---------------------------------------------------------------------------*/
85SYCLInternalLinearAlgebra::ResourceType
86SYCLInternalLinearAlgebra::resource(Matrix const& A)
87{
88 return std::make_tuple(&A.distribution().rowDistribution(),A.blockSize());
89}
90
91void SYCLInternalLinearAlgebra::allocate(ResourceType resource, Vector& v)
92{
93 v.init(*std::get<0>(resource),std::get<1>(resource), true);
94}
95
96void SYCLInternalLinearAlgebra::free(Vector& v)
97{
98 v.clear();
99}
100
101Real SYCLInternalLinearAlgebra::norm0(const SYCLVector<Real>& vx ALIEN_UNUSED_PARAM) const
102{
103 // return SYCLBLASKernel::nrm0(x.space().structInfo(),vx);
104 throw NotImplementedException(
105 A_FUNCINFO, "SYCLLinearAlgebra::norm0 not implemented");
106}
107
108/*---------------------------------------------------------------------------*/
109
110Real SYCLInternalLinearAlgebra::norm1(const SYCLVector<Real>& vx ALIEN_UNUSED_PARAM) const
111{
112 // return SYCLBLASKernel::nrm1(x.space().structInfo(),vx);
113 throw NotImplementedException(
114 A_FUNCINFO, "SYCLLinearAlgebra::norm1 not implemented");
115}
116/*---------------------------------------------------------------------------*/
117
118Real SYCLInternalLinearAlgebra::norm2(const SYCLVector<Real>& vx) const
119{
120#ifdef ALIEN_USE_PERF_TIMER
121 SentryType s(m_timer, "SYCL-NORM2");
122#endif
123 return std::sqrt(m_internal->dot(vx.internal()->values(), vx.internal()->values()));
124}
125
126Real SYCLInternalLinearAlgebra::normInf([[maybe_unused]] const SYCLVector<Real>& vx) const
127{
128#ifdef ALIEN_USE_PERF_TIMER
129 SentryType s(m_timer, "SYCL-NORMINF");
130#endif
131 throw NotImplementedException(
132 A_FUNCINFO, "SYCLLinearAlgebra::normInf not implemented");
133}
134/*---------------------------------------------------------------------------*/
135
136void SYCLInternalLinearAlgebra::mult(const SYCLBEllPackMatrix<Real>& ma,
137 const SYCLVector<Real>& vx,
138 SYCLVector<Real>& vr) const
139{
140#ifdef ALIEN_USE_PERF_TIMER
141 SentryType s(m_timer, "SYCL-SPMV");
142#endif
143 Internal::SYCLBEllPackMatrixMultT<Real>(ma).mult(vx, vr);
144}
145
146void SYCLInternalLinearAlgebra::addLMult(Real alpha,
147 const SYCLBEllPackMatrix<Real>& ma,
148 const SYCLVector<Real>& vx,
149 SYCLVector<Real>& vr) const
150{
151#ifdef ALIEN_USE_PERF_TIMER
152 SentryType s(m_timer, "SYCL-ADDLMULT");
153#endif
154 Internal::SYCLBEllPackMatrixMultT<Real>(ma).addLMult(alpha, vx, vr);
155}
156
157void SYCLInternalLinearAlgebra::addUMult(Real alpha,
158 const SYCLBEllPackMatrix<Real>& ma,
159 const SYCLVector<Real>& vx,
160 SYCLVector<Real>& vr) const
161{
162#ifdef ALIEN_USE_PERF_TIMER
163 SentryType s(m_timer, "SYCL-ADDUMULT");
164#endif
165 Internal::SYCLBEllPackMatrixMultT<Real>(ma).addUMult(alpha, vx, vr);
166}
167
168void SYCLInternalLinearAlgebra::
169multDiag(const SYCLBEllPackMatrix<Real>& ma,
170 const SYCLVector<Real>& vy,
171 SYCLVector<Real>& vr) const
172{
173#ifdef ALIEN_USE_PERF_TIMER
174 SentryType s(m_timer, "SYCL-MULTDIAG");
175#endif
176 Internal::SYCLBEllPackMatrixMultT<Real>(ma).multDiag(vy,vr);
177}
178
179void SYCLInternalLinearAlgebra::
180multDiag(const SYCLVector<Real>& diag,
181 const SYCLVector<Real>& vy,
182 SYCLVector<Real>& vr) const
183{
184#ifdef ALIEN_USE_PERF_TIMER
185 SentryType s(m_timer, "SYCL-MULTDIAG");
186#endif
187 diag.pointWiseMult(vy,vr);
188}
189
190
191void SYCLInternalLinearAlgebra::
192multInvDiag(const SYCLBEllPackMatrix<Real>& ma,
193 SYCLVector<Real>& vr) const
194{
195#ifdef ALIEN_USE_PERF_TIMER
196 SentryType s(m_timer, "SYCL-MULTINVDIAG");
197#endif
198 Internal::SYCLBEllPackMatrixMultT<Real>(ma).multInvDiag(vr);
199}
200
201void SYCLInternalLinearAlgebra::
202computeInvDiag(const SYCLBEllPackMatrix<Real>& ma,
203 SYCLVector<Real>& vr) const
204{
205#ifdef ALIEN_USE_PERF_TIMER
206 SentryType s(m_timer, "SYCL-INVDIAG");
207#endif
208 Internal::SYCLBEllPackMatrixMultT<Real>(ma).computeInvDiag(vr);
209}
210
211/*---------------------------------------------------------------------------*/
212
213void SYCLInternalLinearAlgebra::axpy(Real alpha, const SYCLVector<Real>& vx, SYCLVector<Real>& vr) const
214{
215#ifdef ALIEN_USE_PERF_TIMER
216 SentryType s(m_timer, "SYCL-AXPY");
217#endif
218 m_internal->axpy(alpha, vx.internal()->values(), vr.internal()->values());
219}
220
221/*---------------------------------------------------------------------------*/
222
223void SYCLInternalLinearAlgebra::aypx(Real alpha ALIEN_UNUSED_PARAM,
224 SYCLVector<Real>& y ALIEN_UNUSED_PARAM,
225 const SYCLVector<Real>& x ALIEN_UNUSED_PARAM) const
226{
227 throw NotImplementedException(
228 A_FUNCINFO, "SYCLLinearAlgebra::aypx not implemented");
229}
230
231/*---------------------------------------------------------------------------*/
232
233void SYCLInternalLinearAlgebra::copy(const SYCLVector<Real>& vx, SYCLVector<Real>& vr) const
234{
235#ifdef ALIEN_USE_PERF_TIMER
236 SentryType s(m_timer, "SYCL-COPY");
237#endif
238 m_internal->copy(vx.internal()->values(), vr.internal()->values());
239}
240
241
242void SYCLInternalLinearAlgebra::axpy(Real alpha,
243 const SYCLVector<Real>& vx,
244 Integer stride_x,
245 SYCLVector<Real>& vr,
246 Integer stride_r) const
247{
248#ifdef ALIEN_USE_PERF_TIMER
249 SentryType s(m_timer, "SYCL-AXPY");
250#endif
251 m_internal->axpy(alpha, vx.internal()->values(), stride_x, vr.internal()->values(), stride_r);
252}
253
254/*---------------------------------------------------------------------------*/
255
256void SYCLInternalLinearAlgebra::aypx(Real alpha ALIEN_UNUSED_PARAM,
257 SYCLVector<Real>& y ALIEN_UNUSED_PARAM,
258 Integer stride_y ALIEN_UNUSED_PARAM,
259 const SYCLVector<Real>& x ALIEN_UNUSED_PARAM,
260 Integer stride_x ALIEN_UNUSED_PARAM) const
261{
262 throw NotImplementedException(
263 A_FUNCINFO, "SYCLLinearAlgebra::aypx not implemented");
264}
265
266/*---------------------------------------------------------------------------*/
267
268void SYCLInternalLinearAlgebra::copy(const SYCLVector<Real>& vx,
269 Integer stride_x,
270 SYCLVector<Real>& vr,
271 Integer stride_r) const
272{
273#ifdef ALIEN_USE_PERF_TIMER
274 SentryType s(m_timer, "SYCL-COPY");
275#endif
276 m_internal->copy(vx.internal()->values(), stride_x, vr.internal()->values(), stride_r);
277}
278
279/*---------------------------------------------------------------------------*/
280
281Real SYCLInternalLinearAlgebra::dot(const SYCLVector<Real>& vx, const SYCLVector<Real>& vy) const
282{
283#ifdef ALIEN_USE_PERF_TIMER
284 SentryType s(m_timer, "SYCL-DOT");
285#endif
286 auto value = m_internal->dot(vx.internal()->values(), vy.internal()->values());
287 auto& dist = vx.distribution();
288 if (dist.isParallel()) {
289 return Arccore::MessagePassing::mpAllReduce(dist.parallelMng(),
290 Arccore::MessagePassing::ReduceSum,
291 value);
292 }
293 return value;
294}
295
296void SYCLInternalLinearAlgebra::dot(const SYCLVector<Real>& vx,
297 const SYCLVector<Real>& vy,
298 SYCLInternal::Future<Real>& res) const
299{
300#ifdef ALIEN_USE_PERF_TIMER
301 SentryType s(m_timer, "SYCL-DOT-F");
302#endif
303 m_internal->dot(vx.internal()->values(), vy.internal()->values(), res);
304
305 auto& dist = vx.distribution();
306 if (dist.isParallel()) {
307 res.get() ;
308 Real local_value = res() ;
309 Real* x = &res();
310 auto request = mpNonBlockingAllReduce(dist.parallelMng(),
311 Arccore::MessagePassing::ReduceSum,
312 Arccore::ConstArrayView<Real>(1,&local_value),
313 Arccore::ArrayView<Real>(1,x)) ;
314 res.addRequest(dist.parallelMng(), request);
315 }
316}
317
318/*---------------------------------------------------------------------------*/
319
320void SYCLInternalLinearAlgebra::scal(Real alpha, SYCLVector<Real>& vx) const
321{
322#ifdef ALIEN_USE_PERF_TIMER
323 SentryType s(m_timer, "SYCL-SCAL");
324#endif
325 return m_internal->scal(alpha, vx.internal()->values());
326}
327
328
329void SYCLInternalLinearAlgebra::scal(SYCLVector<Real> const& vx, SYCLBEllPackMatrix<Real>& ma) const
330{
331#ifdef ALIEN_USE_PERF_TIMER
332 SentryType s(m_timer, "SYCL-SCAL");
333#endif
334 return ma.scal(vx);
335}
336
337void SYCLInternalLinearAlgebra::pointwiseMult(const SYCLVector<Real>& vx,
338 const SYCLVector<Real>& vy,
339 SYCLVector<Real>& vz) const
340{
341#ifdef ALIEN_USE_PERF_TIMER
342 SentryType s(m_timer, "SYCL-XYZ");
343#endif
344 return m_internal->pointwiseMult(vx.internal()->values(), vy.internal()->values(), vz.internal()->values());
345}
346
347void SYCLInternalLinearAlgebra::assign(SYCLVector<Real>& vx, Real alpha) const
348{
349#ifdef ALIEN_USE_PERF_TIMER
350 SentryType s(m_timer, "SYCL-ASSIGN");
351#endif
352 return m_internal->assign(alpha, vx.internal()->values());
353}
354
355void SYCLInternalLinearAlgebra::diagonal(
356const SYCLBEllPackMatrix<Real>& ma, SYCLVector<Real>& vx) const
357{
358#ifdef ALIEN_USE_PERF_TIMER
359 SentryType s(m_timer, "SYCL-DIAG");
360#endif
361 Internal::SYCLBEllPackMatrixMultT<Real>(ma).computeDiag(vx);
362}
363
364void SYCLInternalLinearAlgebra::reciprocal(SYCLVector<Real>& x ALIEN_UNUSED_PARAM) const
365{
366 throw NotImplementedException(
367 A_FUNCINFO, "SYCLLinearAlgebra::reciprocal not implemented");
368}
369
370
371Integer SYCLInternalLinearAlgebra::computeCxr(const Matrix& a, Matrix& cxr_a) const
372{
373 auto block_size = a.blockSize() ;
374 cxr_a.setBlockSize(1) ;
375 cxr_a.copy(a) ;
376 return block_size ;
377}
378
379Integer SYCLInternalLinearAlgebra::computeCxr(const Matrix& a, const Vector& diag, Matrix& cxr_a) const
380{
381 auto block_size = a.blockSize() ;
382 cxr_a.setBlockSize(1) ;
383 cxr_a.copy(a) ;
384 cxr_a.scal(diag);
385 return block_size ;
386}
387
388
389
390
391SYCLInternalLinearAlgebraExpr::SYCLInternalLinearAlgebraExpr()
392: IInternalLinearAlgebraExpr<SYCLBEllPackMatrix<Real>, SYCLVector<Real>>()
393{
394 m_internal.reset(new SYCLInternal::KernelInternal());
395}
396
397/*---------------------------------------------------------------------------*/
398
399SYCLInternalLinearAlgebraExpr::~SYCLInternalLinearAlgebraExpr() {}
400
401/*---------------------------------------------------------------------------*/
402/*---------------------------------------------------------------------------*/
403
404Real SYCLInternalLinearAlgebraExpr::norm0(const SYCLVector<Real>& vx ALIEN_UNUSED_PARAM) const
405{
406 // return SYCLBLASKernel::nrm0(x.space().structInfo(),vx);
407 throw NotImplementedException(
408 A_FUNCINFO, "SYCLLinearAlgebra::norm0 not implemented");
409}
410
411/*---------------------------------------------------------------------------*/
412
413Real SYCLInternalLinearAlgebraExpr::norm1(const SYCLVector<Real>& vx ALIEN_UNUSED_PARAM) const
414{
415 // return SYCLBLASKernel::nrm1(x.space().structInfo(),vx);
416 throw NotImplementedException(
417 A_FUNCINFO, "SYCLLinearAlgebra::norm1 not implemented");
418}
419
420/*---------------------------------------------------------------------------*/
421
422Real SYCLInternalLinearAlgebraExpr::norm2(const SYCLVector<Real>& vx) const
423{
424 return std::sqrt(m_internal->dot(vx.internal()->values(), vx.internal()->values()));
425}
426
427Real SYCLInternalLinearAlgebraExpr::normInf([[maybe_unused]] const SYCLVector<Real>& vx) const
428{
429#ifdef ALIEN_USE_PERF_TIMER
430 SentryType s(m_timer, "SYCL-NORMINF");
431#endif
432 throw NotImplementedException(
433 A_FUNCINFO, "SYCLLinearAlgebra::normInf not implemented");
434}
435/*---------------------------------------------------------------------------*/
436
437void SYCLInternalLinearAlgebraExpr::mult(
438const SYCLBEllPackMatrix<Real>& ma, const SYCLVector<Real>& vx, SYCLVector<Real>& vr) const
439{
440 Internal::SYCLBEllPackMatrixMultT<Real>(ma).mult(vx, vr);
441}
442
443void SYCLInternalLinearAlgebraExpr::mult(
444const SYCLBEllPackMatrix<Real>& ma, const UniqueArray<Real>& vx, UniqueArray<Real>& vr) const
445{
446 Internal::SYCLBEllPackMatrixMultT<Real>(ma).mult(vx, vr);
447}
448
449/*---------------------------------------------------------------------------*/
450
451void SYCLInternalLinearAlgebraExpr::axpy([[maybe_unused]] Real alpha,
452 [[maybe_unused]] const UniqueArray<Real>& vx,
453 [[maybe_unused]] UniqueArray<Real>& vr) const
454{
455 throw NotImplementedException(
456 A_FUNCINFO, "SYCLLinearAlgebra::axpy not implemented");
457}
458
459void SYCLInternalLinearAlgebraExpr::axpy(Real alpha, const SYCLVector<Real>& vx, SYCLVector<Real>& vy) const
460{
461 m_internal->axpy(alpha, vx.internal()->values(), vy.internal()->values());
462}
463
464/*---------------------------------------------------------------------------*/
465
466void SYCLInternalLinearAlgebraExpr::aypx([[maybe_unused]] Real alpha, [[maybe_unused]] UniqueArray<Real>& vy, [[maybe_unused]] UniqueArray<Real> const& vx) const
467{
468 throw NotImplementedException(
469 A_FUNCINFO, "SYCLLinearAlgebra::aypx not implemented");
470}
471
472void SYCLInternalLinearAlgebraExpr::aypx(Real alpha ALIEN_UNUSED_PARAM,
473 SYCLVector<Real>& y ALIEN_UNUSED_PARAM, const SYCLVector<Real>& x ALIEN_UNUSED_PARAM) const
474{
475 throw NotImplementedException(
476 A_FUNCINFO, "SYCLLinearAlgebra::aypx not implemented");
477}
478
479/*---------------------------------------------------------------------------*/
480
481void SYCLInternalLinearAlgebraExpr::copy([[maybe_unused]] const UniqueArray<Real>& vx,
482 [[maybe_unused]] UniqueArray<Real>& vr) const
483{
484 throw NotImplementedException(
485 A_FUNCINFO, "SYCLLinearAlgebra::copy not implemented");
486}
487
488void SYCLInternalLinearAlgebraExpr::copy(const SYCLVector<Real>& vx, SYCLVector<Real>& vy) const
489{
490 m_internal->copy(vx.internal()->values(), vy.internal()->values());
491}
492
493/*---------------------------------------------------------------------------*/
494
495Real SYCLInternalLinearAlgebraExpr::dot([[maybe_unused]]
496Integer local_size, [[maybe_unused]] const UniqueArray<Real>& vx, [[maybe_unused]] const UniqueArray<Real>& vy) const
497{
498 throw NotImplementedException(
499 A_FUNCINFO, "SYCLLinearAlgebra::dot not implemented");
500}
501
502Real SYCLInternalLinearAlgebraExpr::dot(const SYCLVector<Real>& vx, const SYCLVector<Real>& vy) const
503{
504 return m_internal->dot(vx.internal()->values(), vy.internal()->values());
505}
506
507Real SYCLInternalLinearAlgebraExpr::norm2([[maybe_unused]] const SYCLBEllPackMatrix<Real>& a) const
508{
509 throw NotImplementedException(
510 A_FUNCINFO, "SYCLLinearAlgebra::notm2 not implemented");
511}
512
513void SYCLInternalLinearAlgebraExpr::copy([[maybe_unused]] const SYCLBEllPackMatrix<Real>& a, [[maybe_unused]] SYCLBEllPackMatrix<Real>& r) const
514{
515 throw NotImplementedException(
516 A_FUNCINFO, "SYCLLinearAlgebra::copy not implemented");
517}
518
519void SYCLInternalLinearAlgebraExpr::add([[maybe_unused]] const SYCLBEllPackMatrix<Real>& a, [[maybe_unused]] SYCLBEllPackMatrix<Real>& r) const
520{
521 throw NotImplementedException(
522 A_FUNCINFO, "SYCLLinearAlgebra::add not implemented");
523}
524
525void SYCLInternalLinearAlgebraExpr::scal([[maybe_unused]] Real alpha, [[maybe_unused]] SYCLBEllPackMatrix<Real>& a) const
526{
527 throw NotImplementedException(
528 A_FUNCINFO, "SYCLLinearAlgebra::scal not implemented");
529}
530/*---------------------------------------------------------------------------*/
531void SYCLInternalLinearAlgebraExpr::scal(Real alpha ALIEN_UNUSED_PARAM, UniqueArray<Real>& x ALIEN_UNUSED_PARAM) const
532{
533 throw NotImplementedException(
534 A_FUNCINFO, "SYCLLinearAlgebra::scal not implemented");
535}
536
537void SYCLInternalLinearAlgebraExpr::scal(Real alpha ALIEN_UNUSED_PARAM, SYCLVector<Real>& x ALIEN_UNUSED_PARAM) const
538{
539 throw NotImplementedException(
540 A_FUNCINFO, "SYCLLinearAlgebra::aypx not implemented");
541}
542
543void SYCLInternalLinearAlgebraExpr::diagonal(
544const SYCLBEllPackMatrix<Real>& a ALIEN_UNUSED_PARAM, SYCLVector<Real>& x ALIEN_UNUSED_PARAM) const
545{
546 throw NotImplementedException(
547 A_FUNCINFO, "SYCLLinearAlgebra::aypx not implemented");
548}
549
550void SYCLInternalLinearAlgebraExpr::reciprocal(SYCLVector<Real>& x ALIEN_UNUSED_PARAM) const
551{
552 throw NotImplementedException(
553 A_FUNCINFO, "SYCLLinearAlgebra::aypx not implemented");
554}
555
556void SYCLInternalLinearAlgebraExpr::pointwiseMult(const SYCLVector<Real>& x ALIEN_UNUSED_PARAM,
557 const SYCLVector<Real>& y ALIEN_UNUSED_PARAM, SYCLVector<Real>& w ALIEN_UNUSED_PARAM) const
558{
559 throw NotImplementedException(
560 A_FUNCINFO, "SYCLLinearAlgebra::aypx not implemented");
561}
562
563template class ALIEN_EXPORT SYCLInternal::Future<double>;
564
565namespace SYCLInternal
566{
567template <>
568sycl::buffer<double>& KernelInternal::getWorkBuffer(std::size_t size)
569{
570 if (m_double_work == nullptr) {
571 m_double_work = new sycl::buffer<double>(size);
572 m_double_work->set_final_data(nullptr);
573 }
574 else {
575 if (size > m_double_work->size()) {
576 delete m_double_work;
577 m_double_work = new sycl::buffer<double>(size);
578 m_double_work->set_final_data(nullptr);
579 }
580 }
581 return *m_double_work;
582}
583}
584
585} // namespace Alien
586
587/*---------------------------------------------------------------------------*/
588/*---------------------------------------------------------------------------*/
LinearAlgebraExprT.h.
LinearAlgebraT.h.
Internal linear algebra interface.
Internal linear algebra interface.
Linear algebra interface.
Linear algebra interface.
-- tab-width: 2; indent-tabs-mode: nil; coding: utf-8-with-signature --
Definition BackEnd.h:17