b2api
B2000++ API Reference Manual, VERSION 4.6
 
Loading...
Searching...
No Matches
b2linalg_solver_slave.H
1//------------------------------------------------------------------------
2// b2linalg_solver_slave.H --
3//
4// written by Mathias Doreille
5// Harald Klimach <harald.klimach@dlr.de>
6// Neda EbrhaimiPour <neda.ebrahimipour@dlr.de>
7//
8// (c) 2011-2016 SMR Engineering & Development SA
9// 2502 Bienne, Switzerland
10//
11// (c) 2021 Deutsches Zentrum für Luft- und Raumfahrt (DLR) e.V.
12// Linder Höhe, 51147 Köln
13//
14// All Rights Reserved. Proprietary source code. The contents of
15// this file may not be disclosed to third parties, copied or
16// duplicated in any form, in whole or in part, without the prior
17// written permission of SMR or DLR.
18//
19//------------------------------------------------------------------------
20
21#ifndef B2LINALG_SOLVER_SLAVE_H_
22#define B2LINALG_SOLVER_SLAVE_H_
23
24//------------------------------------------------------------------------
25//
26// \file
27// Linear algebra slaves for parallel computation.
28//
29//------------------------------------------------------------------------
30
31#include "b2ppconfig.h"
32#include "utils/b2comm.H"
33#include "utils/b2exception.H"
34#include "utils/b2object.H"
35
36#ifdef HAVE_UPMUMPS_DMUMPS_C_H
37#include "MUMPS/dmumps_c.h"
38#else
39#ifdef HAVE_MUMPS_DMUMPS_C_H
40#include "mumps/dmumps_c.h"
41#else
42#ifdef HAVE_DMUMPS_C_H
43#include "dmumps_c.h"
44#endif // HAVE_DMUMPS_C_H
45#endif // HAVE_MUMPS_DMUMPS_C_H
46#endif // HAVE_UPMUMPS_DMUMPS_C_H
47
48#ifdef HAVE_UPMUMPS_ZMUMPS_C_H
49#include "MUMPS/zmumps_c.h"
50#else
51#ifdef HAVE_MUMPS_ZMUMPS_C_H
52#include "mumps/zmumps_c.h"
53#else
54#ifdef HAVE_ZMUMPS_C_H
55#include "zmumps_c.h"
56#endif // HAVE_ZMUMPS_C_H
57#endif // HAVE_MUMPS_ZMUMPS_C_H
58#endif // HAVE_UPMUMPS_ZMUMPS_C_H
59
60#ifdef HAVE_PASTIX_H
61#include <complex>
62extern "C" {
63#include "pastix.h"
64}
65
66#endif // HAVE_PASTIX_H
67
68#include <time.h>
69
70#include <map>
71#include <memory>
72
73#include "utils/b2logging.H"
74
75namespace b2000 { namespace b2linalg {
76
77class SolverSlave {
78public:
79 inline static std::weak_ptr<Communicator> b2comm;
80
81 static void Start() {
82 int buffer;
83 if (logging::get_logger("linear_algebra.sparse_symmetric_solver.mumps")
84 .is_enabled_for(logging::debug)) {
85 buffer = 1;
86 } else {
87 buffer = 0;
88 }
89 std::shared_ptr<Communicator> solvecomm = b2comm.lock();
90 check_comm_ptr(solvecomm, "SolverSlave::Start()");
91 solvecomm->Bcast(&buffer, 1, 0);
92 }
93
94 static void End() {
95 uint64_t buffer[bsize];
96 buffer[0] = 0;
97 std::shared_ptr<Communicator> solvecomm = b2comm.lock();
98 check_comm_ptr(solvecomm, "SolverSlave::End()");
99 solvecomm->Bcast(buffer, bsize);
100 }
101
102 static int GetFcomm() {
103 std::shared_ptr<Communicator> solvecomm = b2comm.lock();
104 check_comm_ptr(solvecomm, "SolverSlave::GetFcomm()");
105 return solvecomm->GetFcomm();
106 }
107
108 static void MumpsSend(unsigned int action, unsigned int type, void* id, unsigned int sym_job) {
109 uint64_t buffer[bsize] = {1, action, type, reinterpret_cast<uint64_t>(id), sym_job};
110 std::shared_ptr<Communicator> solvecomm = b2comm.lock();
111 check_comm_ptr(solvecomm, "SolverSlave::MumpsSend()");
112 solvecomm->Bcast(buffer, bsize);
113 }
114
115 static void PastixSend(unsigned int type, void* id, unsigned int action, uint64_t n) {
116 uint64_t buffer[bsize] = {2, type, reinterpret_cast<uint64_t>(id), action, n};
117 std::shared_ptr<Communicator> solvecomm = b2comm.lock();
118 check_comm_ptr(solvecomm, "SolverSlave::PastixSend()");
119 solvecomm->Bcast(buffer, bsize);
120 }
121
122 void Run() {
123 std::shared_ptr<Communicator> solvecomm = b2comm.lock();
124 check_comm_ptr(solvecomm, "SolverSlave::Run()");
125 {
126 int buffer;
127 solvecomm->Bcast(&buffer, 1, 0);
128 logger_is_enabled_for_debug = buffer == 1 ? true : false;
129 }
130 for (;;) {
131 uint64_t buffer[bsize];
132 solvecomm->Bcast(buffer, bsize);
133
134 switch (buffer[0]) {
135 case 0: // termination of the slave
136 return;
137 case 1: // mumps call
138 Mumps(buffer + 1);
139#ifdef HAVE_PASTIX_H
140 case 2: // pastix call
141 // TODO: PASTIX parallel implementation currently
142 // broken, interface needs to be updated!
143 // Pastix(buffer + 1);
144 return; // termination of the slave, as not usable
145 // with PASTIX.
146#endif // HAVE_PASTIX_H
147 default:
148 assert(1);
149 break;
150 }
151 }
152 }
153
154private:
155 enum StrucType { d = 0, z = 1 };
156 static constexpr int bsize = 5;
157 bool logger_is_enabled_for_debug;
158 std::map<unsigned long int, DMUMPS_STRUC_C> dmumps_object_map;
159 std::map<unsigned long int, ZMUMPS_STRUC_C> zmumps_object_map;
160
161 // Check if the communication object is valid.
162 static void check_comm_ptr(std::shared_ptr<Communicator> comm, std::string label) {
163 if (!comm) {
164 Exception() << "Comm Error: tried to access invalid communicator in"
165 << " " << label << THROW;
166 }
167 }
168
169 void Mumps(uint64_t* buffer) {
170 switch (buffer[0]) {
171 case 1: // new object: action, type, id, sym
172 switch (buffer[1]) {
173 case d: {
174 DMUMPS_STRUC_C& mumps_struct = dmumps_object_map[buffer[2]];
175 mumps_struct.sym = buffer[3];
176 mumps_struct.par = 1;
177 mumps_struct.job = -1;
178 mumps_struct.comm_fortran = static_cast<MUMPS_INT>(GetFcomm());
179 dmumps_c(&mumps_struct);
180 mumps_struct.irn = nullptr;
181 mumps_struct.jcn = nullptr;
182
183 if (logger_is_enabled_for_debug) {
184 mumps_struct.icntl[0] = 6;
185 mumps_struct.icntl[1] = 6;
186 mumps_struct.icntl[2] = 6;
187 mumps_struct.icntl[3] = 3;
188 mumps_struct.icntl[10] = 1;
189 } else {
190 mumps_struct.icntl[0] = 0;
191 mumps_struct.icntl[1] = 0;
192 mumps_struct.icntl[2] = 0;
193 mumps_struct.icntl[3] = 0;
194 }
195 break;
196 }
197 case z: {
198 ZMUMPS_STRUC_C& mumps_struct = zmumps_object_map[buffer[2]];
199 mumps_struct.sym = buffer[3];
200 mumps_struct.par = 1;
201 mumps_struct.job = -1;
202 mumps_struct.comm_fortran = static_cast<MUMPS_INT>(GetFcomm());
203 zmumps_c(&mumps_struct);
204 mumps_struct.irn = nullptr;
205 mumps_struct.jcn = nullptr;
206
207 if (logger_is_enabled_for_debug) {
208 mumps_struct.icntl[0] = 6;
209 mumps_struct.icntl[1] = 6;
210 mumps_struct.icntl[2] = 6;
211 mumps_struct.icntl[3] = 3;
212 mumps_struct.icntl[10] = 1;
213 } else {
214 mumps_struct.icntl[0] = 0;
215 mumps_struct.icntl[1] = 0;
216 mumps_struct.icntl[2] = 0;
217 mumps_struct.icntl[3] = 0;
218 }
219 break;
220 }
221 default:
222 assert(1);
223 break;
224 }
225 break;
226 case 2: // mumps_c: action, type, id, job
227 switch (buffer[1]) {
228 case d: {
229 DMUMPS_STRUC_C& mumps_struct = dmumps_object_map[buffer[2]];
230 mumps_struct.job = buffer[3];
231 dmumps_c(&mumps_struct);
232 break;
233 }
234 case z: {
235 ZMUMPS_STRUC_C& mumps_struct = zmumps_object_map[buffer[2]];
236 mumps_struct.job = buffer[3];
237 zmumps_c(&mumps_struct);
238 break;
239 }
240 default:
241 assert(1);
242 break;
243 }
244 break;
245 default:
246 assert(1);
247 break;
248 }
249 }
250
251#ifdef HAVE_PASTIX_H
252 struct PastixData {
253 // TODO (neda): Adapt to changed PASTIX Interface?
254 // pastix_data_t* pastix_data;
255
256 std::vector<int> ia;
257 std::vector<int> ja;
258 std::vector<double> a;
259 std::vector<int> perm;
260 std::vector<int> invp;
261
262 // TODO (neda): Adapt to changed PASTIX Interface?
263 // pastix_int_t iparm[IPARM_SIZE];
264 // double dparm[DPARM_SIZE];
265 };
266 std::map<unsigned long int, PastixData> pastix_object_map;
267
268 void Pastix(uint64_t* buffer) {
269 PastixData& pd = pastix_object_map[buffer[1]];
270 std::vector<double> rhs(1);
271 std::shared_ptr<Communicator> solvecomm = b2comm.lock();
272 check_comm_ptr(solvecomm, "SolverSlave::Pastix()");
273 switch (buffer[2]) {
274 case 0:
275 // TODO (neda): Adapt to changed PASTIX Interface?
276 // pd.iparm[IPARM_MODIFY_PARAMETER] = API_NO;
277 pd.ia.resize(buffer[3] + 1);
278 pd.perm.resize(buffer[3]);
279 pd.invp.resize(buffer[3]);
280 solvecomm->Bcast(&pd.ia[0], buffer[3] + 1, 0);
281 pd.ja.resize(pd.ia[buffer[3]]);
282 pd.a.resize(pd.ia[buffer[3]]);
283 solvecomm->Bcast(&pd.ja[0], pd.ia[buffer[3]], 0);
284 solvecomm->Bcast(&pd.a[0], pd.ia[buffer[3]], 0);
285 break;
286 case 1:
287 // TODO (neda): Adapt to changed PASTIX Interface?
288 // pd.iparm[IPARM_START_TASK] = API_TASK_ORDERING;
289 // pd.iparm[IPARM_END_TASK] = API_TASK_NUMFACT;
290 // if (logger_is_enabled_for_debug)
291 // pd.iparm[IPARM_VERBOSE] = API_VERBOSE_CHATTERBOX;
292 // else
293 // pd.iparm[IPARM_VERBOSE] = API_VERBOSE_NOT;
294#ifdef HAVE_LIB_TBB
295 // TODO (neda): Adapt to changed PASTIX Interface?
296 // pd.iparm[IPARM_THREAD_NBR] = tbb_default_scheduler
297 // .default_num_threads();
298#endif // HAVE_LIB_TBB
299 // TODO (neda): Adapt to changed PASTIX Interface?
300 // pd.iparm[IPARM_MATRIX_VERIFICATION] = API_YES;
301 break;
302 case 2:
303 // TODO (neda): Adapt to changed PASTIX Interface?
304 // pd.iparm[IPARM_START_TASK] = API_TASK_NUMFACT;
305 // pd.iparm[IPARM_END_TASK] = API_TASK_SOLVE;
306 solvecomm->Bcast(&pd.a[0], pd.ia[buffer[3]], 0);
307 rhs.resize(buffer[3]);
308 solvecomm->Bcast(&rhs[0], buffer[3], 0);
309 break;
310 case 3:
311 // TODO (neda): Adapt to changed PASTIX Interface?
312 // pd.iparm[IPARM_START_TASK] = API_TASK_SOLVE;
313 // pd.iparm[IPARM_END_TASK] = API_TASK_SOLVE;
314 rhs.resize(buffer[3]);
315 solvecomm->Bcast(&rhs[0], buffer[3], 0);
316 break;
317 case 4:
318 // TODO (neda): Adapt to changed PASTIX Interface?
319 // pd.iparm[IPARM_START_TASK] = API_TASK_CLEAN;
320 // pd.iparm[IPARM_END_TASK] = API_TASK_CLEAN;
321 break;
322 default:
323 assert(1);
324 break;
325 }
326 switch (buffer[0]) {
327 // TODO (neda): these interface calls need to be updated.
328 // case d:
329 // d_pastix(&pd.pastix_data, b2comm->get_comm(), buffer[3],
330 // &pd.ia[0], &pd.ja[0], &pd.a[0], &pd.perm[0], &pd.invp[0],
331 // &rhs[0], 1, pd.iparm, pd.dparm);
332 // break;
333 // case z:
334 // z_pastix(&pd.pastix_data, b2comm->get_comm(), 0, 0, 0, 0, 0, 0, 0,
335 // 1, pd.iparm, pd.dparm);
336 // break;
337 default:
338 assert(1);
339 break;
340 }
341 // TODO (neda): these interface calls need to be updated.
342 // if (pd.iparm[IPARM_ERROR_NUMBER] != NO_ERR)
343 // Exception() << "The Pastix solver return error "
344 // << pd.iparm[IPARM_ERROR_NUMBER] << THROW;
345
346 if (buffer[2] == 4) {
347 // clean
348 pastix_object_map.erase(buffer[1]);
349 }
350 }
351#endif // HAVE_PASTIX_H
352}; // clase SolverSlave
353
354}} // namespace b2000::b2linalg
355
356#endif // B2LINALG_SOLVER_SLAVE_H_
#define THROW
Definition b2exception.H:198
Logger & get_logger(const std::string &logger_name="")
Definition b2logging.H:829
Contains the base classes for implementing Finite Elements.
Definition b2boundary_condition.H:32