Tpetra parallel linear algebra Version of the Day
Loading...
Searching...
No Matches
TpetraExt_MatrixMatrix_OpenMP.hpp
1// @HEADER
2// *****************************************************************************
3// Tpetra: Templated Linear Algebra Services Package
4//
5// Copyright 2008 NTESS and the Tpetra contributors.
6// SPDX-License-Identifier: BSD-3-Clause
7// *****************************************************************************
8// @HEADER
9
10#ifndef TPETRA_MATRIXMATRIX_OPENMP_DEF_HPP
11#define TPETRA_MATRIXMATRIX_OPENMP_DEF_HPP
12
13#ifdef HAVE_TPETRA_INST_OPENMP
14namespace Tpetra {
15namespace MMdetails {
16
17template <>
18struct KokkosKernelsSPGEMMBackend<Tpetra::KokkosCompat::KokkosOpenMPWrapperNode> {
19 static std::string parameter_prefix() { return "openmp"; }
20 static std::string algorithm_label() { return "OpenMP"; }
21
22 template <class MatrixType>
23 static void pre_spgemm(MatrixType& Bmerged) {
24#if KOKKOSKERNELS_VERSION >= 50299
25 (void)Bmerged;
26#elif defined(KOKKOSKERNELS_ENABLE_TPL_MKL)
27 using device_t = typename MatrixType::device_type;
28 if constexpr (std::is_same_v<typename device_t::execution_space, Kokkos::OpenMP>) {
29 if (!KokkosSparse::isCrsGraphSorted(Bmerged.graph.row_map, Bmerged.graph.entries)) {
30 Import_Util::sortCrsMatrix(Bmerged);
31 }
32 }
33#else
34 (void)Bmerged;
35#endif
36 }
37
38 template <class Scalar, class LocalOrdinal, class GlobalOrdinal>
39 static void pre_spgemm(const Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode>& A) {
40#if KOKKOSKERNELS_VERSION >= 50299
41 (void)A;
42#elif defined(KOKKOSKERNELS_ENABLE_TPL_MKL)
43 if (A.getCrsGraph()->isSorted()) {
44 return;
45 }
46 auto localA = A.getLocalMatrixDevice();
47 pre_spgemm(localA);
48#else
49 (void)A;
50#endif
51 }
52};
53
54/*********************************************************************************************************/
55// MMM KernelWrappers for Partial Specialization to OpenMP
56template <class Scalar,
57 class LocalOrdinal,
58 class GlobalOrdinal, class LocalOrdinalViewType>
59struct KernelWrappers<Scalar, LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode, LocalOrdinalViewType> {
60 static inline void mult_A_B_newmatrix_kernel_wrapper(CrsMatrixStruct<Scalar, LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode>& Aview,
61 CrsMatrixStruct<Scalar, LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode>& Bview,
62 const LocalOrdinalViewType& Acol2Brow,
63 const LocalOrdinalViewType& Acol2Irow,
64 const LocalOrdinalViewType& Bcol2Ccol,
65 const LocalOrdinalViewType& Icol2Ccol,
66 CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode>& C,
67 Teuchos::RCP<const Import<LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode> > Cimport,
68 const std::string& label = std::string(),
69 const Teuchos::RCP<Teuchos::ParameterList>& params = Teuchos::null) {
70 std::string myalg("SPGEMM_DEFAULT");
71
72 if (!params.is_null()) {
73 if (params->isParameter("openmp: algorithm"))
74 myalg = params->get("openmp: algorithm", myalg);
75 }
76
77 if (myalg == "LTG") {
78 // Use the LTG kernel if requested
79 ::Tpetra::MatrixMatrix::ExtraKernels::mult_A_B_newmatrix_LowThreadGustavsonKernel(Aview, Bview, Acol2Brow, Acol2Irow, Bcol2Ccol, Icol2Ccol, C, Cimport, label, params);
80 } else {
81 Tpetra::MMdetails::kokkos_kernels_mult_A_B_newmatrix(
82 Aview, Bview, Acol2Brow, Acol2Irow, Bcol2Ccol, Icol2Ccol, C, Cimport, label, params);
83 }
84
85#if 0
86 {
87 Teuchos::ArrayRCP< const size_t > Crowptr;
88 Teuchos::ArrayRCP< const LocalOrdinal > Ccolind;
89 Teuchos::ArrayRCP< const Scalar > Cvalues;
90 C.getAllValues(Crowptr,Ccolind,Cvalues);
91
92 // DEBUG
93 int MyPID = C->getComm()->getRank();
94 printf("[%d] Crowptr = ",MyPID);
95 for(size_t i=0; i<(size_t) Crowptr.size(); i++) {
96 printf("%3d ",(int)Crowptr.getConst()[i]);
97 }
98 printf("\n");
99 printf("[%d] Ccolind = ",MyPID);
100 for(size_t i=0; i<(size_t)Ccolind.size(); i++) {
101 printf("%3d ",(int)Ccolind.getConst()[i]);
102 }
103 printf("\n");
104 fflush(stdout);
105 // END DEBUG
106 }
107#endif
108 }
109
110 static inline void mult_A_B_reuse_kernel_wrapper(CrsMatrixStruct<Scalar, LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode>& Aview,
111 CrsMatrixStruct<Scalar, LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode>& Bview,
112 const LocalOrdinalViewType& Acol2Brow,
113 const LocalOrdinalViewType& Acol2Irow,
114 const LocalOrdinalViewType& Bcol2Ccol,
115 const LocalOrdinalViewType& Icol2Ccol,
116 CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode>& C,
117 Teuchos::RCP<const Import<LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode> > Cimport,
118 const std::string& label = std::string(),
119 const Teuchos::RCP<Teuchos::ParameterList>& params = Teuchos::null) {
120#ifdef HAVE_TPETRA_MMM_TIMINGS
121 std::string prefix_mmm = std::string("TpetraExt ") + label + std::string(": ");
122 using Teuchos::TimeMonitor;
123 Teuchos::RCP<TimeMonitor> MM;
124#endif
125
126 // Lots and lots of typedefs
127 using Teuchos::RCP;
128
129 // Options
130 int team_work_size = 16; // Defaults to 16 as per Deveci 12/7/16 - csiefer
131 std::string myalg("LTG");
132 if (!params.is_null()) {
133 if (params->isParameter("openmp: algorithm"))
134 myalg = params->get("openmp: algorithm", myalg);
135 if (params->isParameter("openmp: team work size"))
136 team_work_size = params->get("openmp: team work size", team_work_size);
137 }
138
139 if (myalg == "LTG") {
140 // Use the LTG kernel if requested
141 ::Tpetra::MatrixMatrix::ExtraKernels::mult_A_B_reuse_LowThreadGustavsonKernel(Aview, Bview, Acol2Brow, Acol2Irow, Bcol2Ccol, Icol2Ccol, C, Cimport, label, params);
142 } else {
143 throw std::runtime_error("Tpetra::MatrixMatrix::MMM reuse unknown kernel");
144 }
145
146#ifdef HAVE_TPETRA_MMM_TIMINGS
147 MM = Teuchos::null;
148 MM = rcp(new TimeMonitor(*TimeMonitor::getNewTimer(prefix_mmm + std::string("MMM Reuse OpenMPESFC"))));
149#endif
150 C.fillComplete(C.getDomainMap(), C.getRangeMap());
151 }
152};
153
154// Jacobi KernelWrappers for Partial Specialization to OpenMP
155template <class Scalar,
156 class LocalOrdinal,
157 class GlobalOrdinal, class LocalOrdinalViewType>
158struct KernelWrappers2<Scalar, LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode, LocalOrdinalViewType> {
159 static inline void jacobi_A_B_newmatrix_kernel_wrapper(typename Teuchos::ScalarTraits<Scalar>::magnitudeType omega,
160 const Vector<Scalar, LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode>& Dinv,
161 CrsMatrixStruct<Scalar, LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode>& Aview,
162 CrsMatrixStruct<Scalar, LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode>& Bview,
163 const LocalOrdinalViewType& Acol2Brow,
164 const LocalOrdinalViewType& Acol2Irow,
165 const LocalOrdinalViewType& Bcol2Ccol,
166 const LocalOrdinalViewType& Icol2Ccol,
167 CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode>& C,
168 Teuchos::RCP<const Import<LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode> > Cimport,
169 const std::string& label = std::string(),
170 const Teuchos::RCP<Teuchos::ParameterList>& params = Teuchos::null) {
171#ifdef HAVE_TPETRA_MMM_TIMINGS
172 std::string prefix_mmm = std::string("TpetraExt ") + label + std::string(": ");
173 using Teuchos::TimeMonitor;
174 Teuchos::RCP<TimeMonitor> MM;
175#endif
176
177 // Node-specific code
178 using Teuchos::RCP;
179
180 // Options
181 int team_work_size = 16; // Defaults to 16 as per Deveci 12/7/16 - csiefer
182 std::string myalg("LTG");
183 if (!params.is_null()) {
184 if (params->isParameter("openmp: jacobi algorithm"))
185 myalg = params->get("openmp: jacobi algorithm", myalg);
186 if (params->isParameter("openmp: team work size"))
187 team_work_size = params->get("openmp: team work size", team_work_size);
188 }
189
190 if (myalg == "LTG") {
191 // Use the LTG kernel if requested
192 ::Tpetra::MatrixMatrix::ExtraKernels::jacobi_A_B_newmatrix_LowThreadGustavsonKernel(omega, Dinv, Aview, Bview, Acol2Brow, Acol2Irow, Bcol2Ccol, Icol2Ccol, C, Cimport, label, params);
193 } else if (myalg == "MSAK") {
194 ::Tpetra::MatrixMatrix::ExtraKernels::jacobi_A_B_newmatrix_MultiplyScaleAddKernel(omega, Dinv, Aview, Bview, Acol2Brow, Acol2Irow, Bcol2Ccol, Icol2Ccol, C, Cimport, label, params);
195 } else if (myalg == "KK") {
196 kokkos_kernels_jacobi_A_B_newmatrix(omega, Dinv, Aview, Bview, Acol2Brow, Acol2Irow, Bcol2Ccol, Icol2Ccol, C, Cimport, label, params);
197 } else {
198 throw std::runtime_error("Tpetra::MatrixMatrix::Jacobi newmatrix unknown kernel");
199 }
200
201#ifdef HAVE_TPETRA_MMM_TIMINGS
202 MM = Teuchos::null;
203 MM = rcp(new TimeMonitor(*TimeMonitor::getNewTimer(prefix_mmm + std::string("Jacobi Newmatrix OpenMPESFC"))));
204#endif
205
206 // Final Fillcomplete
207 RCP<Teuchos::ParameterList> labelList = rcp(new Teuchos::ParameterList);
208 labelList->set("Timer Label", label);
209 if (!params.is_null()) labelList->set("compute global constants", params->get("compute global constants", true));
210
211 // NOTE: MSAK already fillCompletes, so we have to check here
212 if (!C.isFillComplete()) {
213 RCP<const Export<LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode> > dummyExport;
214 C.expertStaticFillComplete(Bview.origMatrix->getDomainMap(), Aview.origMatrix->getRangeMap(), Cimport, dummyExport, labelList);
215 }
216 }
217
218 static inline void jacobi_A_B_reuse_kernel_wrapper(typename Teuchos::ScalarTraits<Scalar>::magnitudeType omega,
219 const Vector<Scalar, LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode>& Dinv,
220 CrsMatrixStruct<Scalar, LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode>& Aview,
221 CrsMatrixStruct<Scalar, LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode>& Bview,
222 const LocalOrdinalViewType& Acol2Brow,
223 const LocalOrdinalViewType& Acol2Irow,
224 const LocalOrdinalViewType& Bcol2Ccol,
225 const LocalOrdinalViewType& Icol2Ccol,
226 CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode>& C,
227 Teuchos::RCP<const Import<LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode> > Cimport,
228 const std::string& label = std::string(),
229 const Teuchos::RCP<Teuchos::ParameterList>& params = Teuchos::null) {
230#ifdef HAVE_TPETRA_MMM_TIMINGS
231 std::string prefix_mmm = std::string("TpetraExt ") + label + std::string(": ");
232 using Teuchos::TimeMonitor;
233 Teuchos::RCP<TimeMonitor> MM;
234#endif
235
236 // Lots and lots of typedefs
237 using Teuchos::RCP;
238
239 // Options
240 int team_work_size = 16; // Defaults to 16 as per Deveci 12/7/16 - csiefer
241 std::string myalg("LTG");
242 if (!params.is_null()) {
243 if (params->isParameter("openmp: jacobi algorithm"))
244 myalg = params->get("openmp: jacobi algorithm", myalg);
245 if (params->isParameter("openmp: team work size"))
246 team_work_size = params->get("openmp: team work size", team_work_size);
247 }
248
249 if (myalg == "LTG") {
250 // Use the LTG kernel if requested
251 ::Tpetra::MatrixMatrix::ExtraKernels::jacobi_A_B_reuse_LowThreadGustavsonKernel(omega, Dinv, Aview, Bview, Acol2Brow, Acol2Irow, Bcol2Ccol, Icol2Ccol, C, Cimport, label, params);
252 } else {
253 throw std::runtime_error("Tpetra::MatrixMatrix::Jacobi reuse unknown kernel");
254 }
255
256#ifdef HAVE_TPETRA_MMM_TIMINGS
257 MM = Teuchos::null;
258 MM = rcp(new TimeMonitor(*TimeMonitor::getNewTimer(prefix_mmm + std::string("Jacobi Reuse OpenMPESFC"))));
259#endif
260 C.fillComplete(C.getDomainMap(), C.getRangeMap());
261 }
262};
263
264// Triple-Product KernelWrappers for Partial Specialization to OpenMP
265template <class Scalar,
266 class LocalOrdinal,
267 class GlobalOrdinal, class LocalOrdinalViewType>
268struct KernelWrappers3<Scalar, LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode, LocalOrdinalViewType> {
269 static inline void mult_R_A_P_newmatrix_kernel_wrapper(CrsMatrixStruct<Scalar, LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode>& Rview,
270 CrsMatrixStruct<Scalar, LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode>& Aview,
271 CrsMatrixStruct<Scalar, LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode>& Pview,
272 const LocalOrdinalViewType& Acol2Prow,
273 const LocalOrdinalViewType& Acol2PIrow,
274 const LocalOrdinalViewType& Pcol2Ccol,
275 const LocalOrdinalViewType& PIcol2Ccol,
276 CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode>& Ac,
277 Teuchos::RCP<const Import<LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode> > Acimport,
278 const std::string& label = std::string(),
279 const Teuchos::RCP<Teuchos::ParameterList>& params = Teuchos::null);
280
281 static inline void mult_R_A_P_reuse_kernel_wrapper(CrsMatrixStruct<Scalar, LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode>& Rview,
282 CrsMatrixStruct<Scalar, LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode>& Aview,
283 CrsMatrixStruct<Scalar, LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode>& Pview,
284 const LocalOrdinalViewType& Acol2Prow,
285 const LocalOrdinalViewType& Acol2PIrow,
286 const LocalOrdinalViewType& Pcol2Ccol,
287 const LocalOrdinalViewType& PIcol2Ccol,
288 CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode>& Ac,
289 Teuchos::RCP<const Import<LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode> > Acimport,
290 const std::string& label = std::string(),
291 const Teuchos::RCP<Teuchos::ParameterList>& params = Teuchos::null);
292
293 static inline void mult_PT_A_P_newmatrix_kernel_wrapper(CrsMatrixStruct<Scalar, LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode>& Aview,
294 CrsMatrixStruct<Scalar, LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode>& Pview,
295 const LocalOrdinalViewType& Acol2Prow,
296 const LocalOrdinalViewType& Acol2PIrow,
297 const LocalOrdinalViewType& Pcol2Ccol,
298 const LocalOrdinalViewType& PIcol2Ccol,
299 CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode>& Ac,
300 Teuchos::RCP<const Import<LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode> > Acimport,
301 const std::string& label = std::string(),
302 const Teuchos::RCP<Teuchos::ParameterList>& params = Teuchos::null);
303
304 static inline void mult_PT_A_P_reuse_kernel_wrapper(CrsMatrixStruct<Scalar, LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode>& Aview,
305 CrsMatrixStruct<Scalar, LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode>& Pview,
306 const LocalOrdinalViewType& Acol2Prow,
307 const LocalOrdinalViewType& Acol2PIrow,
308 const LocalOrdinalViewType& Pcol2Ccol,
309 const LocalOrdinalViewType& PIcol2Ccol,
310 CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode>& Ac,
311 Teuchos::RCP<const Import<LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode> > Acimport,
312 const std::string& label = std::string(),
313 const Teuchos::RCP<Teuchos::ParameterList>& params = Teuchos::null);
314};
315
316/*********************************************************************************************************/
317template <class Scalar,
318 class LocalOrdinal,
319 class GlobalOrdinal,
320 class LocalOrdinalViewType>
321void KernelWrappers3<Scalar, LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode, LocalOrdinalViewType>::mult_R_A_P_newmatrix_kernel_wrapper(CrsMatrixStruct<Scalar, LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode>& Rview,
322 CrsMatrixStruct<Scalar, LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode>& Aview,
323 CrsMatrixStruct<Scalar, LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode>& Pview,
324 const LocalOrdinalViewType& Acol2Prow,
325 const LocalOrdinalViewType& Acol2PIrow,
326 const LocalOrdinalViewType& Pcol2Accol,
327 const LocalOrdinalViewType& PIcol2Accol,
328 CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode>& Ac,
329 Teuchos::RCP<const Import<LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode> > Acimport,
330 const std::string& label,
331 const Teuchos::RCP<Teuchos::ParameterList>& params) {
332#ifdef HAVE_TPETRA_MMM_TIMINGS
333 std::string prefix_mmm = std::string("TpetraExt ") + label + std::string(": ");
334 using Teuchos::TimeMonitor;
335 Teuchos::RCP<TimeMonitor> MM;
336#endif
337
338 // Node-specific code
339 std::string nodename("OpenMP");
340
341 // Options
342 std::string myalg("LTG");
343
344 if (!params.is_null()) {
345 if (params->isParameter("openmp: rap algorithm"))
346 myalg = params->get("openmp: rap algorithm", myalg);
347 }
348
349 if (myalg == "LTG") {
350 // Use the LTG kernel if requested
351 ::Tpetra::MatrixMatrix::ExtraKernels::mult_R_A_P_newmatrix_LowThreadGustavsonKernel(Rview, Aview, Pview, Acol2Prow, Acol2PIrow, Pcol2Accol, PIcol2Accol, Ac, Acimport, label, params);
352 } else {
353 throw std::runtime_error("Tpetra::MatrixMatrix::R_A_P newmatrix unknown kernel");
354 }
355}
356
357/*********************************************************************************************************/
358template <class Scalar,
359 class LocalOrdinal,
360 class GlobalOrdinal,
361 class LocalOrdinalViewType>
362void KernelWrappers3<Scalar, LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode, LocalOrdinalViewType>::mult_R_A_P_reuse_kernel_wrapper(
363 CrsMatrixStruct<Scalar, LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode>& Rview,
364 CrsMatrixStruct<Scalar, LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode>& Aview,
365 CrsMatrixStruct<Scalar, LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode>& Pview,
366
367 const LocalOrdinalViewType& Acol2Prow,
368 const LocalOrdinalViewType& Acol2Irow,
369 const LocalOrdinalViewType& Pcol2Ccol,
370 const LocalOrdinalViewType& Icol2Ccol,
371 CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode>& C,
372 Teuchos::RCP<const Import<LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode> > Cimport,
373 const std::string& label,
374 const Teuchos::RCP<Teuchos::ParameterList>& params) {
375#ifdef HAVE_TPETRA_MMM_TIMINGS
376 std::string prefix_mmm = std::string("TpetraExt ") + label + std::string(": ");
377 using Teuchos::TimeMonitor;
378 Teuchos::RCP<TimeMonitor> MM;
379#endif
380
381 // Lots and lots of typedefs
382 using Teuchos::RCP;
383
384 // Options
385 std::string myalg("LTG");
386 if (!params.is_null()) {
387 if (params->isParameter("openmp: rap algorithm"))
388 myalg = params->get("openmp: rap algorithm", myalg);
389 }
390
391 if (myalg == "LTG") {
392 // Use the LTG kernel if requested
393 ::Tpetra::MatrixMatrix::ExtraKernels::mult_R_A_P_reuse_LowThreadGustavsonKernel(Rview, Aview, Pview, Acol2Prow, Acol2Irow, Pcol2Ccol, Icol2Ccol, C, Cimport, label, params);
394 } else {
395 throw std::runtime_error("Tpetra::MatrixMatrix::R_A_P newmatrix unknown kernel");
396 }
397
398#ifdef HAVE_TPETRA_MMM_TIMINGS
399 MM = Teuchos::null;
400 MM = rcp(new TimeMonitor(*TimeMonitor::getNewTimer(prefix_mmm + std::string("RAP Reuse OpenMPESFC"))));
401#endif
402 C.fillComplete(C.getDomainMap(), C.getRangeMap());
403}
404
405/*********************************************************************************************************/
406template <class Scalar,
407 class LocalOrdinal,
408 class GlobalOrdinal,
409 class LocalOrdinalViewType>
410void KernelWrappers3<Scalar, LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode, LocalOrdinalViewType>::mult_PT_A_P_newmatrix_kernel_wrapper(CrsMatrixStruct<Scalar, LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode>& Aview,
411
412 CrsMatrixStruct<Scalar, LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode>& Pview,
413 const LocalOrdinalViewType& Acol2Prow,
414 const LocalOrdinalViewType& Acol2PIrow,
415 const LocalOrdinalViewType& Pcol2Accol,
416 const LocalOrdinalViewType& PIcol2Accol,
417 CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode>& Ac,
418 Teuchos::RCP<const Import<LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode> > Acimport,
419 const std::string& label,
420 const Teuchos::RCP<Teuchos::ParameterList>& params) {
421#ifdef HAVE_TPETRA_MMM_TIMINGS
422 std::string prefix_mmm = std::string("TpetraExt ") + label + std::string(": ");
423 using Teuchos::TimeMonitor;
424 Teuchos::RCP<TimeMonitor> MM;
425#endif
426
427 // Node-specific code
428 std::string nodename("OpenMP");
429
430 // Options
431 std::string myalg("LTG");
432
433 if (!params.is_null()) {
434 if (params->isParameter("openmp: ptap algorithm"))
435 myalg = params->get("openmp: ptap algorithm", myalg);
436 }
437
438 if (myalg == "LTG") {
439#ifdef HAVE_TPETRA_MMM_TIMINGS
440 MM = Teuchos::null;
441 MM = rcp(new TimeMonitor(*TimeMonitor::getNewTimer(prefix_mmm + std::string("PTAP local transpose"))));
442#endif
443
444 using Teuchos::ParameterList;
445 using Teuchos::RCP;
446 using LO = LocalOrdinal;
447 using GO = GlobalOrdinal;
448 using SC = Scalar;
449
450 // We don't need a kernel-level PTAP, we just transpose here
451 using transposer_type =
452 RowMatrixTransposer<SC, LO, GO, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode>;
453 transposer_type transposer(Pview.origMatrix, label + "XP: ");
454 RCP<ParameterList> transposeParams(new ParameterList);
455 if (!params.is_null()) {
456 transposeParams->set("compute global constants",
457 params->get("compute global constants: temporaries",
458 false));
459 }
460 transposeParams->set("sort", false);
461 RCP<CrsMatrix<SC, LO, GO, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode> > Ptrans =
462 transposer.createTransposeLocal(transposeParams);
463 CrsMatrixStruct<SC, LO, GO, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode> Rview;
464 Rview.origMatrix = Ptrans;
465
466 using ::Tpetra::MatrixMatrix::ExtraKernels::mult_R_A_P_newmatrix_LowThreadGustavsonKernel;
467 mult_R_A_P_newmatrix_LowThreadGustavsonKernel(Rview, Aview, Pview, Acol2Prow, Acol2PIrow, Pcol2Accol,
468 PIcol2Accol, Ac, Acimport, label, params);
469 } else {
470 throw std::runtime_error("Tpetra::MatrixMatrix::PT_A_P newmatrix unknown kernel");
471 }
472}
473
474/*********************************************************************************************************/
475template <class Scalar,
476 class LocalOrdinal,
477 class GlobalOrdinal,
478 class LocalOrdinalViewType>
479void KernelWrappers3<Scalar, LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode, LocalOrdinalViewType>::mult_PT_A_P_reuse_kernel_wrapper(CrsMatrixStruct<Scalar, LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode>& Aview,
480
481 CrsMatrixStruct<Scalar, LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode>& Pview,
482 const LocalOrdinalViewType& Acol2Prow,
483 const LocalOrdinalViewType& Acol2PIrow,
484 const LocalOrdinalViewType& Pcol2Accol,
485 const LocalOrdinalViewType& PIcol2Accol,
486 CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode>& Ac,
487 Teuchos::RCP<const Import<LocalOrdinal, GlobalOrdinal, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode> > Acimport,
488 const std::string& label,
489 const Teuchos::RCP<Teuchos::ParameterList>& params) {
490#ifdef HAVE_TPETRA_MMM_TIMINGS
491 std::string prefix_mmm = std::string("TpetraExt ") + label + std::string(": ");
492 using Teuchos::TimeMonitor;
493 Teuchos::RCP<TimeMonitor> MM;
494#endif
495
496 // Node-specific code
497 std::string nodename("OpenMP");
498
499 // Options
500 std::string myalg("LTG");
501
502 if (!params.is_null()) {
503 if (params->isParameter("openmp: ptap algorithm"))
504 myalg = params->get("openmp: ptap algorithm", myalg);
505 }
506
507 if (myalg == "LTG") {
508#ifdef HAVE_TPETRA_MMM_TIMINGS
509 MM = Teuchos::null;
510 MM = rcp(new TimeMonitor(*TimeMonitor::getNewTimer(prefix_mmm + std::string("PTAP local transpose"))));
511#endif
512
513 using Teuchos::ParameterList;
514 using Teuchos::RCP;
515 using LO = LocalOrdinal;
516 using GO = GlobalOrdinal;
517 using SC = Scalar;
518
519 // We don't need a kernel-level PTAP, we just transpose here
520 using transposer_type =
521 RowMatrixTransposer<SC, LO, GO, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode>;
522 transposer_type transposer(Pview.origMatrix, label + "XP: ");
523 RCP<ParameterList> transposeParams(new ParameterList);
524 if (!params.is_null()) {
525 transposeParams->set("compute global constants",
526 params->get("compute global constants: temporaries",
527 false));
528 }
529 transposeParams->set("sort", false);
530 RCP<CrsMatrix<SC, LO, GO, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode> > Ptrans =
531 transposer.createTransposeLocal(transposeParams);
532 CrsMatrixStruct<SC, LO, GO, Tpetra::KokkosCompat::KokkosOpenMPWrapperNode> Rview;
533 Rview.origMatrix = Ptrans;
534
535 using ::Tpetra::MatrixMatrix::ExtraKernels::mult_R_A_P_reuse_LowThreadGustavsonKernel;
536 mult_R_A_P_reuse_LowThreadGustavsonKernel(Rview, Aview, Pview, Acol2Prow, Acol2PIrow, Pcol2Accol,
537 PIcol2Accol, Ac, Acimport, label, params);
538 } else {
539 throw std::runtime_error("Tpetra::MatrixMatrix::PT_A_P reuse unknown kernel");
540 }
541 Ac.fillComplete(Ac.getDomainMap(), Ac.getRangeMap());
542}
543
544} // namespace MMdetails
545} // namespace Tpetra
546
547#endif // OpenMP
548
549#endif
Tpetra::CrsMatrix::getLocalMatrixDevice
TPETRA_DETAILS_ALWAYS_INLINE local_matrix_device_type getLocalMatrixDevice() const
The local sparse matrix.
Definition Tpetra_CrsMatrix_def.hpp:1011
Tpetra::CrsMatrix::getCrsGraph
Teuchos::RCP< const crs_graph_type > getCrsGraph() const
This matrix's graph, as a CrsGraph.
Definition Tpetra_CrsMatrix_def.hpp:978
Tpetra
Namespace Tpetra contains the class and methods constituting the Tpetra library.
Generated by doxygen 1.15.0