Alien  1.3.0
Developer documentation
Loading...
Searching...
No Matches
SimpleCSRInternalLinearAlgebra.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 <cmath>
10#include <tuple>
11
12#include "SimpleCSRInternalLinearAlgebra.h"
13#include <alien/utils/Precomp.h>
14
15#include <arccore/base/NotImplementedException.h>
16#include <arccore/base/TraceInfo.h>
17
18#include <alien/core/backend/LinearAlgebraExprT.h>
19#include <alien/core/backend/LinearAlgebraT.h>
20
21#include "CBLASMPIKernel.h"
22#include "SimpleCSRMatrixMult.h"
23
24/*---------------------------------------------------------------------------*/
25
26namespace Alien
27{
28
29using namespace Arccore;
30
31template class ALIEN_EXPORT LinearAlgebra<BackEnd::tag::simplecsr>;
32template class ALIEN_EXPORT LinearAlgebraExpr<BackEnd::tag::simplecsr>;
33
34ALIEN_EXPORT IInternalLinearAlgebra<SimpleCSRMatrix<Real>, SimpleCSRVector<Real>>*
35SimpleCSRInternalLinearAlgebraFactory()
36{
37 return new SimpleCSRInternalLinearAlgebra();
38}
39
40ALIEN_EXPORT IInternalLinearAlgebraExpr<SimpleCSRMatrix<Real>, SimpleCSRVector<Real>>*
41SimpleCSRInternalLinearAlgebraExprFactory()
42{
43 return new SimpleCSRInternalLinearAlgebraExpr();
44}
45
46/*---------------------------------------------------------------------------*/
47
48namespace Internal = SimpleCSRInternal;
49
50/*---------------------------------------------------------------------------*/
51
52SimpleCSRInternalLinearAlgebra::SimpleCSRInternalLinearAlgebra()
53: IInternalLinearAlgebra<CSRMatrix, CSRVector>()
54{}
55
56/*---------------------------------------------------------------------------*/
57
58SimpleCSRInternalLinearAlgebra::~SimpleCSRInternalLinearAlgebra()
59{
60#ifdef ALIEN_USE_PERF_TIMER
61 m_timer.printInfo("SIMPLECSR-ALGEBRA");
62#endif
63}
64
65/*---------------------------------------------------------------------------*/
66/*---------------------------------------------------------------------------*/
67SimpleCSRInternalLinearAlgebra::ResourceType
68SimpleCSRInternalLinearAlgebra::resource(Matrix const& A)
69{
70 return std::make_tuple(&A.distribution().rowDistribution(),A.blockSize());
71}
72
73void SimpleCSRInternalLinearAlgebra::allocate(ResourceType resource, Vector& v)
74{
75 v.init(*std::get<0>(resource),std::get<1>(resource), true);
76}
77
78void SimpleCSRInternalLinearAlgebra::free(Vector& v)
79{
80 v.clear();
81}
82
83Real SimpleCSRInternalLinearAlgebra::norm0(const CSRVector& vx) const
84{
85#ifdef ALIEN_USE_PERF_TIMER
86 SentryType s(m_timer, "CSR-NORM0");
87#endif
88
89 return CBLASMPIKernel::nrm0(vx.distribution(), vx);
90}
91
92/*---------------------------------------------------------------------------*/
93
94Real SimpleCSRInternalLinearAlgebra::norm1(const CSRVector& vx) const
95{
96#ifdef ALIEN_USE_PERF_TIMER
97 SentryType s(m_timer, "CSR-NORM1");
98#endif
99
100 return CBLASMPIKernel::nrm1(vx.distribution(), vx);
101}
102
103/*---------------------------------------------------------------------------*/
104
105Real SimpleCSRInternalLinearAlgebra::norm2(const CSRVector& vx) const
106{
107#ifdef ALIEN_USE_PERF_TIMER
108 SentryType s(m_timer, "CSR-NORM2");
109#endif
110 return CBLASMPIKernel::nrm2(vx.distribution(), vx);
111}
112
113/*---------------------------------------------------------------------------*/
114
115Real SimpleCSRInternalLinearAlgebra::normInf(const CSRVector& vx) const
116{
117#ifdef ALIEN_USE_PERF_TIMER
118 SentryType s(m_timer, "CSR-NORMINF");
119#endif
120 return CBLASMPIKernel::nrmInf(vx.distribution(), vx);
121}
122
123/*---------------------------------------------------------------------------*/
124
125void SimpleCSRInternalLinearAlgebra::synchronize(const CSRMatrix& ma,
126 CSRVector& vx) const
127{
128 Internal::SimpleCSRMatrixMultT<Real>(ma).synchronize(vx);
129}
130
131void SimpleCSRInternalLinearAlgebra::mult(const CSRMatrix& ma,
132 const CSRVector& vx,
133 CSRVector& vr) const
134{
135#ifdef ALIEN_USE_PERF_TIMER
136 SentryType s(m_timer, "CSR-SPMV");
137#endif
138 Internal::SimpleCSRMatrixMultT<Real>(ma).mult(vx, vr);
139}
140
141void SimpleCSRInternalLinearAlgebra::addLMult(Real alpha,
142 const CSRMatrix& ma,
143 const CSRVector& vx,
144 CSRVector& vr) const
145{
146#ifdef ALIEN_USE_PERF_TIMER
147 SentryType s(m_timer, "CSR-AddLMult");
148#endif
149 Internal::SimpleCSRMatrixMultT<Real>(ma).addLMult(alpha, vx, vr);
150}
151
152void SimpleCSRInternalLinearAlgebra::addUMult(Real alpha,
153 const CSRMatrix& ma,
154 const CSRVector& vx,
155 CSRVector& vr) const
156{
157#ifdef ALIEN_USE_PERF_TIMER
158 SentryType s(m_timer, "CSR-AddUMult");
159#endif
160 Internal::SimpleCSRMatrixMultT<Real>(ma).addUMult(alpha, vx, vr);
161}
162
163void SimpleCSRInternalLinearAlgebra::multDiag(const CSRMatrix& ma, CSRVector& vr) const
164{
165#ifdef ALIEN_USE_PERF_TIMER
166 SentryType s(m_timer, "CSR-MULTDIAG");
167#endif
168 Internal::SimpleCSRMatrixMultT<Real>(ma).multDiag(vr);
169}
170
171void SimpleCSRInternalLinearAlgebra::computeDiag(const CSRMatrix& ma, CSRVector& vr) const
172{
173#ifdef ALIEN_USE_PERF_TIMER
174 SentryType s(m_timer, "CSR-DIAG");
175#endif
176 Internal::SimpleCSRMatrixMultT<Real>(ma).computeDiag(vr);
177}
178
179void SimpleCSRInternalLinearAlgebra::multInvDiag(const CSRMatrix& ma, CSRVector& vr) const
180{
181#ifdef ALIEN_USE_PERF_TIMER
182 SentryType s(m_timer, "CSR-MULTINVDIAG");
183#endif
184 Internal::SimpleCSRMatrixMultT<Real>(ma).multInvDiag(vr);
185}
186
187void SimpleCSRInternalLinearAlgebra::computeInvDiag(const CSRMatrix& ma, CSRVector& vr) const
188{
189#ifdef ALIEN_USE_PERF_TIMER
190 SentryType s(m_timer, "CSR-INVDIAG");
191#endif
192 Internal::SimpleCSRMatrixMultT<Real>(ma).computeInvDiag(vr);
193}
194
195/*---------------------------------------------------------------------------*/
196
197void SimpleCSRInternalLinearAlgebra::axpy(Real alpha, const CSRVector& vx, CSRVector& vr) const
198{
199#ifdef ALIEN_USE_PERF_TIMER
200 SentryType s(m_timer, "CSR-AXPY");
201#endif
202 CBLASMPIKernel::axpy(vx.distribution(), alpha, vx, vr);
203}
204
205/*---------------------------------------------------------------------------*/
206
207void SimpleCSRInternalLinearAlgebra::aypx(Real alpha ALIEN_UNUSED_PARAM,
208 CSRVector& y ALIEN_UNUSED_PARAM,
209 const CSRVector& x ALIEN_UNUSED_PARAM) const
210{
211 throw NotImplementedException(
212 A_FUNCINFO, "SimpleCSRLinearAlgebra::aypx not implemented");
213}
214
215/*---------------------------------------------------------------------------*/
216
217void SimpleCSRInternalLinearAlgebra::copy(const CSRVector& vx, CSRVector& vr) const
218{
219#ifdef ALIEN_USE_PERF_TIMER
220 SentryType s(m_timer, "CSR-COPY");
221#endif
222 CBLASMPIKernel::copy(vx.distribution(), vx, vr);
223}
224
225void SimpleCSRInternalLinearAlgebra::axpy(Real alpha,
226 const CSRVector& vx,
227 Integer stride_x,
228 CSRVector& vr,
229 Integer stride_r) const
230{
231#ifdef ALIEN_USE_PERF_TIMER
232 SentryType s(m_timer, "CSR-AXPY");
233#endif
234 CBLASMPIKernel::axpy(vx.distribution(), alpha, vx, stride_x, vr, stride_r);
235}
236
237/*---------------------------------------------------------------------------*/
238
239void SimpleCSRInternalLinearAlgebra::aypx([[maybe_unused]] Real alpha,
240 [[maybe_unused]] CSRVector& y,
241 [[maybe_unused]] Integer stride_y,
242 [[maybe_unused]] const CSRVector& x,
243 [[maybe_unused]] Integer stride_x) const
244{
245 throw NotImplementedException(
246 A_FUNCINFO, "SimpleCSRLinearAlgebra::aypx not implemented");
247}
248
249void SimpleCSRInternalLinearAlgebra::copy(const Vector& vx,
250 Integer stride_x,
251 Vector& vr,
252 Integer stride_r) const
253{
254 CBLASMPIKernel::copy(vx.distribution(), vx, stride_x, vr, stride_r);
255}
256/*---------------------------------------------------------------------------*/
257
258Real SimpleCSRInternalLinearAlgebra::dot(const CSRVector& vx, const CSRVector& vy) const
259{
260#ifdef ALIEN_USE_PERF_TIMER
261 SentryType s(m_timer, "CSR-DOT");
262#endif
263 return CBLASMPIKernel::dot(vx.distribution(), vx, vy);
264}
265
266void SimpleCSRInternalLinearAlgebra::dot(const CSRVector& vx,
267 const CSRVector& vy,
268 SimpleCSRInternalLinearAlgebra::FutureType& res) const
269{
270#ifdef ALIEN_USE_PERF_TIMER
271 SentryType s(m_timer, "CSR-DOT-F");
272#endif
273 res() = CBLASMPIKernel::dot(vx.distribution(), vx, vy);
274}
275
276/*---------------------------------------------------------------------------*/
277
278void SimpleCSRInternalLinearAlgebra::scal(Real alpha, CSRVector& vx) const
279{
280#ifdef ALIEN_USE_PERF_TIMER
281 SentryType s(m_timer, "CSR-SCAL");
282#endif
283 CBLASMPIKernel::scal(vx.distribution(), alpha, vx);
284}
285
286void SimpleCSRInternalLinearAlgebra::scal(CSRVector const& vx, CSRMatrix& a) const
287{
288#ifdef ALIEN_USE_PERF_TIMER
289 SentryType s(m_timer, "CSR-SCAL");
290#endif
291 a.scal(vx.data()) ;
292}
293
294void SimpleCSRInternalLinearAlgebra::diagonal(
295const CSRMatrix& a, CSRVector& diag) const
296{
297#ifdef ALIEN_USE_PERF_TIMER
298 SentryType s(m_timer, "CSR-DIAG");
299#endif
300 Internal::SimpleCSRMatrixMultT<Real>(a).computeDiag(diag);
301}
302
303void SimpleCSRInternalLinearAlgebra::reciprocal(CSRVector& x ALIEN_UNUSED_PARAM) const
304{
305 throw NotImplementedException(
306 A_FUNCINFO, "SimpleCSRLinearAlgebra::reciprocal not implemented");
307}
308
309void SimpleCSRInternalLinearAlgebra::pointwiseMult(const CSRVector& vx,
310 const CSRVector& vy,
311 CSRVector& vz) const
312{
313#ifdef ALIEN_USE_PERF_TIMER
314 SentryType s(m_timer, "CSR-XYZ");
315#endif
316 CBLASMPIKernel::pointwiseMult(vx.distribution(), vx, vy, vz);
317}
318
319void SimpleCSRInternalLinearAlgebra::assign(CSRVector& vx, Real alpha) const
320{
321#ifdef ALIEN_USE_PERF_TIMER
322 SentryType s(m_timer, "CSR-ASSIGN");
323#endif
324 CBLASMPIKernel::assign(vx.distribution(), alpha, vx);
325}
326
327Integer SimpleCSRInternalLinearAlgebra::computeCxr(const CSRMatrix& a, CSRMatrix& cxr_a) const
328{
329 auto block_size = a.blockSize() ;
330 cxr_a.setBlockSize(1) ;
331 cxr_a.copy(a) ;
332 return block_size ;
333}
334
335Integer SimpleCSRInternalLinearAlgebra::computeCxr(const CSRMatrix& a, const CSRVector& diag, CSRMatrix& cxr_a) const
336{
337 auto block_size = a.blockSize() ;
338 cxr_a.setBlockSize(1) ;
339 cxr_a.copy(a) ;
340 cxr_a.scal(diag.data());
341 return block_size ;
342}
343
344
345SimpleCSRInternalLinearAlgebraExpr::SimpleCSRInternalLinearAlgebraExpr()
346: IInternalLinearAlgebraExpr<CSRMatrix, CSRVector>()
347{}
348
349/*---------------------------------------------------------------------------*/
350
351SimpleCSRInternalLinearAlgebraExpr::~SimpleCSRInternalLinearAlgebraExpr() {}
352
353/*---------------------------------------------------------------------------*/
354/*---------------------------------------------------------------------------*/
355
356Real SimpleCSRInternalLinearAlgebraExpr::norm0(const CSRVector& vx) const
357{
358 return CBLASMPIKernel::nrm0(vx.distribution(), vx);
359}
360
361/*---------------------------------------------------------------------------*/
362
363Real SimpleCSRInternalLinearAlgebraExpr::norm1(const CSRVector& vx) const
364{
365 return CBLASMPIKernel::nrm1(vx.distribution(), vx);
366}
367
368/*---------------------------------------------------------------------------*/
369
370Real SimpleCSRInternalLinearAlgebraExpr::norm2(const CSRVector& vx) const
371{
372 return CBLASMPIKernel::nrm2(vx.distribution(), vx);
373}
374
375/*---------------------------------------------------------------------------*/
376
377Real SimpleCSRInternalLinearAlgebraExpr::normInf(const CSRVector& vx) const
378{
379 return CBLASMPIKernel::nrmInf(vx.distribution(), vx);
380}
381
382/*---------------------------------------------------------------------------*/
383
384Real SimpleCSRInternalLinearAlgebraExpr::norm2(const CSRMatrix& mx) const
385{
386#ifdef ALIEN_USE_PERF_TIMER
387 SentryType s(m_timer, "MATRIX-CSR-NORM2");
388#endif
389 return CBLASMPIKernel::matrix_nrm2(mx.distribution(), mx);
390}
391/*---------------------------------------------------------------------------*/
392
393void SimpleCSRInternalLinearAlgebraExpr::mult(
394const CSRMatrix& ma, const CSRVector& vx, CSRVector& vr) const
395{
396 Internal::SimpleCSRMatrixMultT<Real>(ma).mult(vx, vr);
397}
398
399void SimpleCSRInternalLinearAlgebraExpr::mult(
400const CSRMatrix& ma, const UniqueArray<Real>& vx, UniqueArray<Real>& vr) const
401{
402 Internal::SimpleCSRMatrixMultT<Real>(ma).mult(vx, vr);
403}
404
405/*---------------------------------------------------------------------------*/
406
407void SimpleCSRInternalLinearAlgebraExpr::axpy(Real alpha, const UniqueArray<Real>& vx, UniqueArray<Real>& vr) const
408{
409 cblas::axpy(vx.size(), alpha, dataPtr(vx), 1, dataPtr(vr), 1);
410}
411
412void SimpleCSRInternalLinearAlgebraExpr::axpy(Real alpha, const CSRVector& vx, CSRVector& vr) const
413{
414 CBLASMPIKernel::axpy(vx.distribution(), alpha, vx, vr);
415}
416
417/*---------------------------------------------------------------------------*/
418
419void SimpleCSRInternalLinearAlgebraExpr::aypx(Real alpha, UniqueArray<Real>& vy, UniqueArray<Real> const& vx) const
420{
421 throw NotImplementedException(
422 A_FUNCINFO, "SimpleCSRLinearAlgebra::aypx not implemented");
423}
424
425void SimpleCSRInternalLinearAlgebraExpr::aypx(Real alpha ALIEN_UNUSED_PARAM,
426 CSRVector& y ALIEN_UNUSED_PARAM, const CSRVector& x ALIEN_UNUSED_PARAM) const
427{
428 throw NotImplementedException(
429 A_FUNCINFO, "SimpleCSRLinearAlgebra::aypx not implemented");
430}
431
432/*---------------------------------------------------------------------------*/
433
434void SimpleCSRInternalLinearAlgebraExpr::copy(
435const UniqueArray<Real>& vx, UniqueArray<Real>& vr) const
436{
437 cblas::copy(vx.size(), dataPtr(vx), 1, dataPtr(vr), 1);
438}
439
440void SimpleCSRInternalLinearAlgebraExpr::copy(const CSRVector& vx, CSRVector& vr) const
441{
442 CBLASMPIKernel::copy(vx.distribution(), vx, vr);
443}
444
445void SimpleCSRInternalLinearAlgebraExpr::copy(const CSRMatrix& ma, CSRMatrix& mr) const
446{
447#ifdef ALIEN_USE_PERF_TIMER
448 SentryType s(m_timer, "MATRIX-CSR-COPY");
449#endif
450 mr.copy(ma);
451}
452
453void SimpleCSRInternalLinearAlgebraExpr::add(const CSRMatrix& ma, CSRMatrix& mr) const
454{
455#ifdef ALIEN_USE_PERF_TIMER
456 SentryType s(m_timer, "MATRIX-CSR-ADD");
457#endif
458
459 cblas::axpy(mr.getProfile().getNnz(), 1, (CSRMatrix::ValueType*)ma.data(), 1, mr.data(), 1);
460}
461
462void SimpleCSRInternalLinearAlgebraExpr::scal(Real alpha, CSRMatrix& mr) const
463{
464#ifdef ALIEN_USE_PERF_TIMER
465 SentryType s(m_timer, "MATRIX-CSR-SCAL");
466#endif
467
468 cblas::scal(mr.getProfile().getNnz(), alpha, mr.data(), 1);
469}
470
471/*---------------------------------------------------------------------------*/
472
473Real SimpleCSRInternalLinearAlgebraExpr::dot(
474Integer local_size, const UniqueArray<Real>& vx, const UniqueArray<Real>& vy) const
475{
476 return cblas::dot(local_size, dataPtr(vx), 1, dataPtr(vy), 1);
477}
478
479Real SimpleCSRInternalLinearAlgebraExpr::dot(const CSRVector& vx, const CSRVector& vy) const
480{
481 return CBLASMPIKernel::dot(vx.distribution(), vx, vy);
482}
483
484/*---------------------------------------------------------------------------*/
485void SimpleCSRInternalLinearAlgebraExpr::scal(Real alpha ALIEN_UNUSED_PARAM, UniqueArray<Real>& x ALIEN_UNUSED_PARAM) const
486{
487 throw NotImplementedException(
488 A_FUNCINFO, "SimpleCSRLinearAlgebra::scal not implemented");
489}
490
491void SimpleCSRInternalLinearAlgebraExpr::scal(Real alpha ALIEN_UNUSED_PARAM, CSRVector& x ALIEN_UNUSED_PARAM) const
492{
493 throw NotImplementedException(
494 A_FUNCINFO, "SimpleCSRLinearAlgebra::aypx not implemented");
495}
496
497void SimpleCSRInternalLinearAlgebraExpr::diagonal(
498const CSRMatrix& a ALIEN_UNUSED_PARAM, CSRVector& x ALIEN_UNUSED_PARAM) const
499{
500 throw NotImplementedException(
501 A_FUNCINFO, "SimpleCSRLinearAlgebra::aypx not implemented");
502}
503
504void SimpleCSRInternalLinearAlgebraExpr::reciprocal(CSRVector& x ALIEN_UNUSED_PARAM) const
505{
506 throw NotImplementedException(
507 A_FUNCINFO, "SimpleCSRLinearAlgebra::aypx not implemented");
508}
509
510void SimpleCSRInternalLinearAlgebraExpr::pointwiseMult(const CSRVector& x ALIEN_UNUSED_PARAM,
511 const CSRVector& y ALIEN_UNUSED_PARAM, CSRVector& w ALIEN_UNUSED_PARAM) const
512{
513 throw NotImplementedException(
514 A_FUNCINFO, "SimpleCSRLinearAlgebra::aypx not implemented");
515}
516
517} // namespace Alien
518
519/*---------------------------------------------------------------------------*/
520/*---------------------------------------------------------------------------*/
LinearAlgebraExprT.h
LinearAlgebraExprT.h.
LinearAlgebraT.h
LinearAlgebraT.h.
Alien::IInternalLinearAlgebraExpr
Internal linear algebra interface.
Definition IInternalLinearAlgebraExprT.h:39
Alien::IInternalLinearAlgebra
Internal linear algebra interface.
Definition IInternalLinearAlgebraT.h:34
Alien::LinearAlgebraExpr
Linear algebra interface.
Definition LinearAlgebraExpr.h:49
Alien::LinearAlgebra
Linear algebra interface.
Definition LinearAlgebra.h:50
Alien
-- tab-width: 2; indent-tabs-mode: nil; coding: utf-8-with-signature --
Definition BackEnd.h:17