Alien  1.3.0
User 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
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
35SimpleCSRInternalLinearAlgebraFactory()
36{
38}
39
41SimpleCSRInternalLinearAlgebraExprFactory()
42{
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::lowerUpperSynchronize(const CSRMatrix& ma,
164 CSRVector& vr) const
165{
166#ifdef ALIEN_USE_PERF_TIMER
167 SentryType s(m_timer, "CSR-luSynchronize");
168#endif
169 Internal::SimpleCSRMatrixMultT<Real>(ma).lowerUpperSynchronize(vr, vr);
170}
171
172void SimpleCSRInternalLinearAlgebra::upperLowerSynchronize(const CSRMatrix& ma,
173 CSRVector& vr) const
174{
175#ifdef ALIEN_USE_PERF_TIMER
176 SentryType s(m_timer, "CSR-ulSynchronize");
177#endif
178 Internal::SimpleCSRMatrixMultT<Real>(ma).upperLowerSynchronize(vr, vr);
179}
180
181void SimpleCSRInternalLinearAlgebra::multDiag(const CSRMatrix& ma, CSRVector& vr) const
182{
183#ifdef ALIEN_USE_PERF_TIMER
184 SentryType s(m_timer, "CSR-MULTDIAG");
185#endif
186 Internal::SimpleCSRMatrixMultT<Real>(ma).multDiag(vr);
187}
188
189void SimpleCSRInternalLinearAlgebra::computeDiag(const CSRMatrix& ma, CSRVector& vr) const
190{
191#ifdef ALIEN_USE_PERF_TIMER
192 SentryType s(m_timer, "CSR-DIAG");
193#endif
194 Internal::SimpleCSRMatrixMultT<Real>(ma).computeDiag(vr);
195}
196
197void SimpleCSRInternalLinearAlgebra::multInvDiag(const CSRMatrix& ma, CSRVector& vr) const
198{
199#ifdef ALIEN_USE_PERF_TIMER
200 SentryType s(m_timer, "CSR-MULTINVDIAG");
201#endif
202 Internal::SimpleCSRMatrixMultT<Real>(ma).multInvDiag(vr);
203}
204
205void SimpleCSRInternalLinearAlgebra::computeInvDiag(const CSRMatrix& ma, CSRVector& vr) const
206{
207#ifdef ALIEN_USE_PERF_TIMER
208 SentryType s(m_timer, "CSR-INVDIAG");
209#endif
210 Internal::SimpleCSRMatrixMultT<Real>(ma).computeInvDiag(vr);
211}
212
213/*---------------------------------------------------------------------------*/
214
215void SimpleCSRInternalLinearAlgebra::axpy(Real alpha, const CSRVector& vx, CSRVector& vr) const
216{
217#ifdef ALIEN_USE_PERF_TIMER
218 SentryType s(m_timer, "CSR-AXPY");
219#endif
220 CBLASMPIKernel::axpy(vx.distribution(), alpha, vx, vr);
221}
222
223/*---------------------------------------------------------------------------*/
224
225void SimpleCSRInternalLinearAlgebra::aypx(Real alpha ALIEN_UNUSED_PARAM,
226 CSRVector& y ALIEN_UNUSED_PARAM,
227 const CSRVector& x ALIEN_UNUSED_PARAM) const
228{
229 throw NotImplementedException(
230 A_FUNCINFO, "SimpleCSRLinearAlgebra::aypx not implemented");
231}
232
233/*---------------------------------------------------------------------------*/
234
235void SimpleCSRInternalLinearAlgebra::copy(const CSRVector& vx, CSRVector& vr) const
236{
237#ifdef ALIEN_USE_PERF_TIMER
238 SentryType s(m_timer, "CSR-COPY");
239#endif
240 CBLASMPIKernel::copy(vx.distribution(), vx, vr);
241}
242
243void SimpleCSRInternalLinearAlgebra::axpy(Real alpha,
244 const CSRVector& vx,
245 Integer stride_x,
246 CSRVector& vr,
247 Integer stride_r) const
248{
249#ifdef ALIEN_USE_PERF_TIMER
250 SentryType s(m_timer, "CSR-AXPY");
251#endif
252 CBLASMPIKernel::axpy(vx.distribution(), alpha, vx, stride_x, vr, stride_r);
253}
254
255/*---------------------------------------------------------------------------*/
256
257void SimpleCSRInternalLinearAlgebra::aypx([[maybe_unused]] Real alpha,
258 [[maybe_unused]] CSRVector& y,
259 [[maybe_unused]] Integer stride_y,
260 [[maybe_unused]] const CSRVector& x,
261 [[maybe_unused]] Integer stride_x) const
262{
263 throw NotImplementedException(
264 A_FUNCINFO, "SimpleCSRLinearAlgebra::aypx not implemented");
265}
266
267void SimpleCSRInternalLinearAlgebra::copy(const Vector& vx,
268 Integer stride_x,
269 Vector& vr,
270 Integer stride_r) const
271{
272 CBLASMPIKernel::copy(vx.distribution(), vx, stride_x, vr, stride_r);
273}
274/*---------------------------------------------------------------------------*/
275
276Real SimpleCSRInternalLinearAlgebra::dot(const CSRVector& vx, const CSRVector& vy) const
277{
278#ifdef ALIEN_USE_PERF_TIMER
279 SentryType s(m_timer, "CSR-DOT");
280#endif
281 return CBLASMPIKernel::dot(vx.distribution(), vx, vy);
282}
283
284void SimpleCSRInternalLinearAlgebra::dot(const CSRVector& vx,
285 const CSRVector& vy,
286 SimpleCSRInternalLinearAlgebra::FutureType& res) const
287{
288#ifdef ALIEN_USE_PERF_TIMER
289 SentryType s(m_timer, "CSR-DOT-F");
290#endif
291 res() = CBLASMPIKernel::dot(vx.distribution(), vx, vy);
292}
293
294/*---------------------------------------------------------------------------*/
295
296void SimpleCSRInternalLinearAlgebra::scal(Real alpha, CSRVector& vx) const
297{
298#ifdef ALIEN_USE_PERF_TIMER
299 SentryType s(m_timer, "CSR-SCAL");
300#endif
301 CBLASMPIKernel::scal(vx.distribution(), alpha, vx);
302}
303
304void SimpleCSRInternalLinearAlgebra::scal(CSRVector const& vx, CSRMatrix& a) const
305{
306#ifdef ALIEN_USE_PERF_TIMER
307 SentryType s(m_timer, "CSR-SCAL");
308#endif
309 a.scal(vx.data()) ;
310}
311
312void SimpleCSRInternalLinearAlgebra::diagonal(
313const CSRMatrix& a, CSRVector& diag) const
314{
315#ifdef ALIEN_USE_PERF_TIMER
316 SentryType s(m_timer, "CSR-DIAG");
317#endif
318 Internal::SimpleCSRMatrixMultT<Real>(a).computeDiag(diag);
319}
320
321void SimpleCSRInternalLinearAlgebra::reciprocal(CSRVector& x ALIEN_UNUSED_PARAM) const
322{
323 throw NotImplementedException(
324 A_FUNCINFO, "SimpleCSRLinearAlgebra::reciprocal not implemented");
325}
326
327void SimpleCSRInternalLinearAlgebra::pointwiseMult(const CSRVector& vx,
328 const CSRVector& vy,
329 CSRVector& vz) const
330{
331#ifdef ALIEN_USE_PERF_TIMER
332 SentryType s(m_timer, "CSR-XYZ");
333#endif
334 CBLASMPIKernel::pointwiseMult(vx.distribution(), vx, vy, vz);
335}
336
337void SimpleCSRInternalLinearAlgebra::assign(CSRVector& vx, Real alpha) const
338{
339#ifdef ALIEN_USE_PERF_TIMER
340 SentryType s(m_timer, "CSR-ASSIGN");
341#endif
342 CBLASMPIKernel::assign(vx.distribution(), alpha, vx);
343}
344
345Integer SimpleCSRInternalLinearAlgebra::computeCxr(const CSRMatrix& a, CSRMatrix& cxr_a) const
346{
347 auto block_size = a.blockSize() ;
348 cxr_a.setBlockSize(1) ;
349 cxr_a.copy(a) ;
350 return block_size ;
351}
352
353Integer SimpleCSRInternalLinearAlgebra::computeCxr(const CSRMatrix& a, const CSRVector& diag, CSRMatrix& cxr_a) const
354{
355 auto block_size = a.blockSize() ;
356 cxr_a.setBlockSize(1) ;
357 cxr_a.copy(a) ;
358 cxr_a.scal(diag.data());
359 return block_size ;
360}
361
362
363SimpleCSRInternalLinearAlgebraExpr::SimpleCSRInternalLinearAlgebraExpr()
364: IInternalLinearAlgebraExpr<CSRMatrix, CSRVector>()
365{}
366
367/*---------------------------------------------------------------------------*/
368
369SimpleCSRInternalLinearAlgebraExpr::~SimpleCSRInternalLinearAlgebraExpr() {}
370
371/*---------------------------------------------------------------------------*/
372/*---------------------------------------------------------------------------*/
373
374Real SimpleCSRInternalLinearAlgebraExpr::norm0(const CSRVector& vx) const
375{
376 return CBLASMPIKernel::nrm0(vx.distribution(), vx);
377}
378
379/*---------------------------------------------------------------------------*/
380
381Real SimpleCSRInternalLinearAlgebraExpr::norm1(const CSRVector& vx) const
382{
383 return CBLASMPIKernel::nrm1(vx.distribution(), vx);
384}
385
386/*---------------------------------------------------------------------------*/
387
388Real SimpleCSRInternalLinearAlgebraExpr::norm2(const CSRVector& vx) const
389{
390 return CBLASMPIKernel::nrm2(vx.distribution(), vx);
391}
392
393/*---------------------------------------------------------------------------*/
394
395Real SimpleCSRInternalLinearAlgebraExpr::normInf(const CSRVector& vx) const
396{
397 return CBLASMPIKernel::nrmInf(vx.distribution(), vx);
398}
399
400/*---------------------------------------------------------------------------*/
401
402Real SimpleCSRInternalLinearAlgebraExpr::norm2(const CSRMatrix& mx) const
403{
404#ifdef ALIEN_USE_PERF_TIMER
405 SentryType s(m_timer, "MATRIX-CSR-NORM2");
406#endif
407 return CBLASMPIKernel::matrix_nrm2(mx.distribution(), mx);
408}
409/*---------------------------------------------------------------------------*/
410
411void SimpleCSRInternalLinearAlgebraExpr::mult(
412const CSRMatrix& ma, const CSRVector& vx, CSRVector& vr) const
413{
414 Internal::SimpleCSRMatrixMultT<Real>(ma).mult(vx, vr);
415}
416
417void SimpleCSRInternalLinearAlgebraExpr::mult(
418const CSRMatrix& ma, const UniqueArray<Real>& vx, UniqueArray<Real>& vr) const
419{
420 Internal::SimpleCSRMatrixMultT<Real>(ma).mult(vx, vr);
421}
422
423/*---------------------------------------------------------------------------*/
424
425void SimpleCSRInternalLinearAlgebraExpr::axpy(Real alpha, const UniqueArray<Real>& vx, UniqueArray<Real>& vr) const
426{
427 cblas::axpy(vx.size(), alpha, dataPtr(vx), 1, dataPtr(vr), 1);
428}
429
430void SimpleCSRInternalLinearAlgebraExpr::axpy(Real alpha, const CSRVector& vx, CSRVector& vr) const
431{
432 CBLASMPIKernel::axpy(vx.distribution(), alpha, vx, vr);
433}
434
435/*---------------------------------------------------------------------------*/
436
437void SimpleCSRInternalLinearAlgebraExpr::aypx(Real alpha, UniqueArray<Real>& vy, UniqueArray<Real> const& vx) const
438{
439 throw NotImplementedException(
440 A_FUNCINFO, "SimpleCSRLinearAlgebra::aypx not implemented");
441}
442
443void SimpleCSRInternalLinearAlgebraExpr::aypx(Real alpha ALIEN_UNUSED_PARAM,
444 CSRVector& y ALIEN_UNUSED_PARAM, const CSRVector& x ALIEN_UNUSED_PARAM) const
445{
446 throw NotImplementedException(
447 A_FUNCINFO, "SimpleCSRLinearAlgebra::aypx not implemented");
448}
449
450/*---------------------------------------------------------------------------*/
451
452void SimpleCSRInternalLinearAlgebraExpr::copy(
453const UniqueArray<Real>& vx, UniqueArray<Real>& vr) const
454{
455 cblas::copy(vx.size(), dataPtr(vx), 1, dataPtr(vr), 1);
456}
457
458void SimpleCSRInternalLinearAlgebraExpr::copy(const CSRVector& vx, CSRVector& vr) const
459{
460 CBLASMPIKernel::copy(vx.distribution(), vx, vr);
461}
462
463void SimpleCSRInternalLinearAlgebraExpr::copy(const CSRMatrix& ma, CSRMatrix& mr) const
464{
465#ifdef ALIEN_USE_PERF_TIMER
466 SentryType s(m_timer, "MATRIX-CSR-COPY");
467#endif
468 mr.copy(ma);
469}
470
471void SimpleCSRInternalLinearAlgebraExpr::add(const CSRMatrix& ma, CSRMatrix& mr) const
472{
473#ifdef ALIEN_USE_PERF_TIMER
474 SentryType s(m_timer, "MATRIX-CSR-ADD");
475#endif
476
477 cblas::axpy(mr.getProfile().getNnz(), 1, (CSRMatrix::ValueType*)ma.data(), 1, mr.data(), 1);
478}
479
480void SimpleCSRInternalLinearAlgebraExpr::scal(Real alpha, CSRMatrix& mr) const
481{
482#ifdef ALIEN_USE_PERF_TIMER
483 SentryType s(m_timer, "MATRIX-CSR-SCAL");
484#endif
485
486 cblas::scal(mr.getProfile().getNnz(), alpha, mr.data(), 1);
487}
488
489/*---------------------------------------------------------------------------*/
490
491Real SimpleCSRInternalLinearAlgebraExpr::dot(
492Integer local_size, const UniqueArray<Real>& vx, const UniqueArray<Real>& vy) const
493{
494 return cblas::dot(local_size, dataPtr(vx), 1, dataPtr(vy), 1);
495}
496
497Real SimpleCSRInternalLinearAlgebraExpr::dot(const CSRVector& vx, const CSRVector& vy) const
498{
499 return CBLASMPIKernel::dot(vx.distribution(), vx, vy);
500}
501
502/*---------------------------------------------------------------------------*/
503void SimpleCSRInternalLinearAlgebraExpr::scal(Real alpha ALIEN_UNUSED_PARAM, UniqueArray<Real>& x ALIEN_UNUSED_PARAM) const
504{
505 throw NotImplementedException(
506 A_FUNCINFO, "SimpleCSRLinearAlgebra::scal not implemented");
507}
508
509void SimpleCSRInternalLinearAlgebraExpr::scal(Real alpha ALIEN_UNUSED_PARAM, CSRVector& x ALIEN_UNUSED_PARAM) const
510{
511 throw NotImplementedException(
512 A_FUNCINFO, "SimpleCSRLinearAlgebra::aypx not implemented");
513}
514
515void SimpleCSRInternalLinearAlgebraExpr::diagonal(
516const CSRMatrix& a ALIEN_UNUSED_PARAM, CSRVector& x ALIEN_UNUSED_PARAM) const
517{
518 throw NotImplementedException(
519 A_FUNCINFO, "SimpleCSRLinearAlgebra::aypx not implemented");
520}
521
522void SimpleCSRInternalLinearAlgebraExpr::reciprocal(CSRVector& x ALIEN_UNUSED_PARAM) const
523{
524 throw NotImplementedException(
525 A_FUNCINFO, "SimpleCSRLinearAlgebra::aypx not implemented");
526}
527
528void SimpleCSRInternalLinearAlgebraExpr::pointwiseMult(const CSRVector& x ALIEN_UNUSED_PARAM,
529 const CSRVector& y ALIEN_UNUSED_PARAM, CSRVector& w ALIEN_UNUSED_PARAM) const
530{
531 throw NotImplementedException(
532 A_FUNCINFO, "SimpleCSRLinearAlgebra::aypx not implemented");
533}
534
535} // namespace Alien
536
537/*---------------------------------------------------------------------------*/
538/*---------------------------------------------------------------------------*/
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