Alien  1.3.0
User documentation
Loading...
Searching...
No Matches
SYCLBEllPackMatrixMultT.h
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#pragma once
9
10#include <arccore/collections/Array2.h>
11
12#include <alien/kernels/sycl/data/SYCLSendRecvOp.h>
13/*---------------------------------------------------------------------------*/
14
15namespace Alien::SYCLInternal
16{
17
18/*---------------------------------------------------------------------------*/
19
20template <typename ValueT>
21SYCLBEllPackMatrixMultT<ValueT>::SYCLBEllPackMatrixMultT(const MatrixType& matrix)
22: m_matrix_impl(matrix)
23{}
24
25/*---------------------------------------------------------------------------*/
26
27template <typename ValueT>
28void SYCLBEllPackMatrixMultT<ValueT>::mult(const VectorType& x, VectorType& y) const
29{
30 if (m_matrix_impl.m_is_parallel)
31 _parallelMult(x, y);
32 else
33 _seqMult(x, y);
34}
35
36template <typename ValueT>
37void SYCLBEllPackMatrixMultT<ValueT>::addLMult(Real alpha, const VectorType& x, VectorType& y) const
38{
39#ifdef ALIEN_USE_PERF_TIMER
40 typename MatrixType::SentryType sentry(m_matrix_impl.timer(), "SYCL-AddLMult");
41#endif
42
43 //if (m_matrix_impl.m_is_parallel)
44 // lowerUpperSynchronize(x,const_cast<VectorType&>(x)) ;
45
46
47 m_matrix_impl.addLMult(alpha, x, y);
48
49 //if (m_matrix_impl.m_is_parallel)
50 // m_matrix_impl.endDistLMult(alpha, x, y) ;
51}
52
53template <typename ValueT>
54void SYCLBEllPackMatrixMultT<ValueT>::addUMult(Real alpha, const VectorType& x, VectorType& y) const
55{
56#ifdef ALIEN_USE_PERF_TIMER
57 typename MatrixType::SentryType sentry(m_matrix_impl.timer(), "SYCL-AddUMult");
58#endif
59 //if (m_matrix_impl.m_is_parallel)
60 // lowerUpperSynchronize(x,const_cast<VectorType&>(x)) ;
61
62 m_matrix_impl.addUMult(alpha, x, y);
63
64 //if (m_matrix_impl.m_is_parallel)
65 // m_matrix_impl.endDistUMult(alpha, x, y) ;
66}
67
68template <typename ValueT>
69void SYCLBEllPackMatrixMultT<ValueT>::mult(const UniqueArray<Real>& x, UniqueArray<Real>& y) const
70{
71#ifdef ALIEN_USE_PERF_TIMER
72 typename MatrixType::SentryType sentry(m_matrix_impl.timer(), "SYCL-SPMV");
73#endif
74 if (m_matrix_impl.m_is_parallel)
75 _parallelMult(x, y);
76 else
77 _seqMult(x, y);
78}
79
80/*---------------------------------------------------------------------------*/
81
82template <typename ValueT>
83void SYCLBEllPackMatrixMultT<ValueT>::_parallelMult(
84const VectorType& x_impl, VectorType& y_impl) const
85{
86 //Alien::alien_debug([&] {Alien::cout() << "SYCL PARALLEL MULT : "<<m_matrix_impl.getGhostSize();});
87 //Universe().traceMng()->flush() ;
88 int block_size = m_matrix_impl.blockSize() ;
89
90 SYCLSendRecvOp<ValueT> op(x_impl.internal()->values(),
91 m_matrix_impl.m_matrix_dist_info.m_send_info,
92 m_matrix_impl.internal()->getSendIds(),
93 m_matrix_impl.m_send_policy,
94 x_impl.internal()->ghostValues(m_matrix_impl.getGhostSize()*block_size),
95 m_matrix_impl.m_matrix_dist_info.m_recv_info,
96 m_matrix_impl.internal()->getRecvIds(),
97 m_matrix_impl.m_recv_policy,
98 m_matrix_impl.m_parallel_mng,
99 m_matrix_impl.m_trace,
100 block_size);
101
102 op.start();
103
104 m_matrix_impl.mult(x_impl, y_impl);
105
106 op.end();
107
108 m_matrix_impl.endDistMult(x_impl, y_impl);
109
110 //Alien::alien_debug([&] {Alien::cout() << "End SYCL PARALLEL MULT";});
111 //Universe().traceMng()->flush() ;
112}
113
114template <typename ValueT>
115void SYCLBEllPackMatrixMultT<ValueT>::_parallelMult(
116[[maybe_unused]] const UniqueArray<Real>& x_impl, [[maybe_unused]] UniqueArray<Real>& y_impl) const
117{
118#ifdef ENABLE_MPI_SYCL
119 Real* y_ptr = dataPtr(y_impl);
120 Real* x_ptr = (Real*)dataPtr(x_impl);
121 ConstArrayView<Real> matrix = m_matrix_impl.m_matrix.getValues();
122 ConstArrayView<Integer> cols = m_matrix_impl.getDistStructInfo().m_cols;
123 ConstArrayView<Integer> row_offset =
124 m_matrix_impl.m_matrix.getProfile().getRowOffset();
125 SendRecvOp<Real> op(x_ptr, m_matrix_impl.m_matrix_dist_info.m_send_info,
126 m_matrix_impl.m_send_policy, x_ptr, m_matrix_impl.m_matrix_dist_info.m_recv_info,
127 m_matrix_impl.m_recv_policy, m_matrix_impl.m_parallel_mng, m_matrix_impl.m_trace);
128 op.start();
129#endif
130
131 //m_matrix_impl.mult(x_impl,y_impl) ;
132
133#ifdef ENABLE_MPI_SYCL
134 op.end();
135
136 Integer interface_nrow = m_matrix_impl.m_matrix_dist_info.m_interface_nrow;
137 ConstArrayView<Integer> row_ids = m_matrix_impl.m_matrix_dist_info.m_interface_rows;
138 for (Integer i = 0; i < interface_nrow; ++i) {
139 Integer irow = row_ids[i];
140 Integer off = row_offset[irow] + local_row_size[irow];
141 Integer off2 = row_offset[irow + 1];
142 Real tmpy = 0.;
143 for (Integer j = off; j < off2; ++j) {
144 tmpy += matrix[j] * x_ptr[cols[j]];
145 }
146 y_ptr[irow] += tmpy;
147 }
148#endif
149}
150/*---------------------------------------------------------------------------*/
151
152template <typename ValueT>
153void SYCLBEllPackMatrixMultT<ValueT>::_seqMult(const VectorType& x_impl, VectorType& y_impl) const
154{
155 m_matrix_impl.mult(x_impl, y_impl);
156}
157
158template <typename ValueT>
159void SYCLBEllPackMatrixMultT<ValueT>::_seqMult([[maybe_unused]] const UniqueArray<Real>& x_impl,
160 [[maybe_unused]] UniqueArray<Real>& y_impl) const
161{
162}
163
164template <typename ValueT>
165void SYCLBEllPackMatrixMultT<ValueT>::lowerUpperSynchronize(const VectorType& x_impl, VectorType& y_impl) const
166{
167 int block_size = m_matrix_impl.blockSize() ;
168 alien_debug([&]{cout()<<"BEGIN lowerUpperSynchronize "<<block_size;}) ;
169
170 SYCLSendRecvOp<ValueT> op(x_impl.internal()->values(),
171 m_matrix_impl.m_matrix_dist_info.m_send_info,
172 m_matrix_impl.internal()->getSendIds(),
173 m_matrix_impl.m_send_policy,
174 y_impl.internal()->ghostValues(m_matrix_impl.getGhostSize()*block_size),
175 m_matrix_impl.m_matrix_dist_info.m_recv_info,
176 m_matrix_impl.internal()->getRecvIds(),
177 m_matrix_impl.m_recv_policy,
178 m_matrix_impl.m_parallel_mng,
179 m_matrix_impl.m_trace,
180 block_size);
181
182 op.lowerRecv();
183 op.upperSend();
184 alien_debug([&]{cout()<<"AFTER lowerUpperSynchronize "<<block_size;}) ;
185}
186
187template <typename ValueT>
188void SYCLBEllPackMatrixMultT<ValueT>::upperLowerSynchronize(const VectorType& x_impl, VectorType& y_impl) const
189{
190 int block_size = m_matrix_impl.blockSize() ;
191 alien_debug([&]{cout()<<"BEGIN upperLowerSynchronize "<<block_size;}) ;
192
193 SYCLSendRecvOp<ValueT> op(x_impl.internal()->values(),
194 m_matrix_impl.m_matrix_dist_info.m_send_info,
195 m_matrix_impl.internal()->getSendIds(),
196 m_matrix_impl.m_send_policy,
197 y_impl.internal()->ghostValues(m_matrix_impl.getGhostSize()*block_size),
198 m_matrix_impl.m_matrix_dist_info.m_recv_info,
199 m_matrix_impl.internal()->getRecvIds(),
200 m_matrix_impl.m_recv_policy,
201 m_matrix_impl.m_parallel_mng,
202 m_matrix_impl.m_trace,
203 block_size);
204
205 op.upperRecv();
206 op.lowerSend();
207 alien_debug([&]{cout()<<"END upperLowerSynchronize "<<block_size;}) ;
208}
209
210
211/*---------------------------------------------------------------------------*/
212template <typename ValueT>
213void SYCLBEllPackMatrixMultT<ValueT>::multDiag(VectorType const& y, VectorType& z) const
214{
215 m_matrix_impl.multDiag(y,z);
216}
217
218template <typename ValueT>
219void SYCLBEllPackMatrixMultT<ValueT>::multDiag(VectorType& y) const
220{
221 m_matrix_impl.multDiag(y);
222}
223
224template <typename ValueT>
225void SYCLBEllPackMatrixMultT<ValueT>::computeDiag(VectorType& y) const
226{
227 m_matrix_impl.computeDiag(y);
228}
229
230template <typename ValueT>
231void SYCLBEllPackMatrixMultT<ValueT>::multInvDiag(VectorType& y) const
232{
233 m_matrix_impl.multInvDiag(y);
234}
235
236template <typename ValueT>
237void SYCLBEllPackMatrixMultT<ValueT>::computeInvDiag(VectorType& y) const
238{
239 m_matrix_impl.computeInvDiag(y);
240}
241
242/*---------------------------------------------------------------------------*/
243
244/*---------------------------------------------------------------------------*/
245
246} // namespace Alien::SYCLInternal
247
248/*---------------------------------------------------------------------------*/
Alien::SYCLInternal::SYCLBEllPackMatrixMultT::mult
void mult(const VectorType &x, VectorType &y) const
Matrix vector product.
Definition SYCLBEllPackMatrixMultT.h:28
Alien::SYCLInternal::SYCLBEllPackMatrixMultT::SYCLBEllPackMatrixMultT
SYCLBEllPackMatrixMultT(const MatrixType &matrix)
Constructeur de la classe.
Definition SYCLBEllPackMatrixMultT.h:21