PaStiX Handbook  6.2.1
solver.h
Go to the documentation of this file.
1 /**
2  * @file solver.h
3  *
4  * PaStiX solver structure header.
5  *
6  * @copyright 2004-2021 Bordeaux INP, CNRS (LaBRI UMR 5800), Inria,
7  * Univ. Bordeaux. All rights reserved.
8  *
9  * @version 6.2.0
10  * @author David Goudin
11  * @author Pascal Henon
12  * @author Francois Pellegrini
13  * @author Pierre Ramet
14  * @author Mathieu Faverge
15  * @author Xavier Lacoste
16  * @author Esragul Korkmaz
17  * @author Gregoire Pichon
18  * @author Tony Delarue
19  * @date 2021-03-30
20  *
21  **/
22 #ifndef _solver_h_
23 #define _solver_h_
24 
25 struct blendctrl_s;
26 typedef struct blendctrl_s BlendCtrl;
27 
28 struct simuctrl_s;
29 typedef struct simuctrl_s SimuCtrl;
30 
31 #include "pastix_lowrank.h"
32 
33 /**
34  * @name Cblk properties
35  * @{
36  * The type and structure definitions.
37  * Define the mask for the cblks in the cblktype field:
38  * - 1st bit: The cblk is a local fanin accumulating the contributions for a remote cblk
39  * - 2nd bit: The cblk is stored in a 2D layout fashion as in a tiled matrix, otherwise the standard 1D lapack layout is used
40  * - 3rd bit: The cblk generates 2D granularity tasks, instead of a single 1D tasks that perform factorization, solves and updates
41  * - 4th bit: The cblk is compressed in Low-Rank (implies CBLK_LAYOUT_2D), otherwise it is stored in dense
42  * - 5th bit: The cblk is part of the Schur complement if set
43  * - 6th bit: The cblk is a local copy of a remote fanin to be received.
44  */
45 #define CBLK_FANIN (1 << 0)
46 #define CBLK_LAYOUT_2D (1 << 1)
47 #define CBLK_TASKS_2D (1 << 2)
48 #define CBLK_COMPRESSED (1 << 3)
49 #define CBLK_IN_SCHUR (1 << 4)
50 #define CBLK_IN_LAST (1 << 5)
51 #define CBLK_RECV (1 << 6)
52 
53 /**
54  *@}
55  */
56 
57 /*
58  * The type and structure definitions.
59  */
60 #define COMP_1D 0
61 #define DIAG 1
62 #define E1 2
63 #define E2 3
64 #define DRUNK 4
65 
66 /**
67  * @brief Solver recv block structure.
68  */
69 typedef struct solver_blok_recv_s {
70  pastix_int_t frownum; /**< First row contributed in the block (lrownum+1 of the original block if none) */
71  pastix_int_t lrownum; /**< Last row contributed in the block (frownum-1 of the original block if none) */
73 
74 /**
75  * @brief Solver recv column block structure.
76  *
77  * The structure works as a list of cblk for each remote owner to compute the exact contributions.
78  */
79 typedef struct solver_cblk_recv_s {
80  struct solver_cblk_recv_s *next;
81  pastix_int_t ownerid;
82  pastix_int_t fcolnum; /**< First column contributed in the block (lcolnum+1 of the original cblk if none) */
83  pastix_int_t lcolnum; /**< Last column contributed in the block (fcolnum-1 of the original cblk if none) */
84  solver_blok_recv_t bloktab[1]; /**< Array of reception blocks */
86 
87 /**
88  * @brief The task structure for the numerical factorization
89  */
90 typedef struct task_s {
91  pastix_int_t taskid; /**< COMP_1D DIAG E1 E2 */
92  pastix_int_t prionum; /**< Priority value for the factorization */
93  pastix_int_t cblknum; /**< Attached column block */
94  pastix_int_t bloknum; /**< Attached block */
95  pastix_int_t volatile ctrbcnt; /**< Total number of contributions */
96 #if defined(PASTIX_DYNSCHED)
97  int threadid;/**< Index of the bubble which contains the task */
98 #endif
99 } Task;
100 
101 #define GPUID_UNDEFINED -2 /**< GPU still undefined */
102 #define GPUID_NONE -1 /**< Block not computed on GPU */
103 
104 /**
105  * @brief Solver block structure.
106  */
107 typedef struct solver_blok_s {
108  void *handler[2]; /**< Runtime data handler */
109  pastix_int_t lcblknm; /**< Local column block */
110  pastix_int_t fcblknm; /**< Facing column block */
111  pastix_int_t gbloknm; /**< Index in global bloktab (UNUSED) */
112  pastix_int_t frownum; /**< First row index */
113  pastix_int_t lrownum; /**< Last row index (inclusive) */
114  pastix_int_t coefind; /**< Index in coeftab */
115  pastix_int_t browind; /**< Index in browtab */
116  int8_t gpuid; /**< Store on which GPU the block is computed */
117  int8_t inlast; /**< Index of the block among last separator (2), coupling with last separator (1) or other blocks (0) */
118  int iluklvl; /**< The block ILU(k) level */
119 
120  /* LR structures */
121  pastix_lrblock_t *LRblock[2]; /**< Store the blok (L/U) in LR format. Allocated for the cblk. */
122 } SolverBlok;
123 
124 /**
125  * @brief Solver column block structure.
126  */
127 typedef struct solver_cblk_s {
128  pastix_atomic_lock_t lock; /**< Lock to protect computation on the cblk */
129  volatile int32_t ctrbcnt; /**< Number of contribution to receive */
130  int8_t cblktype; /**< Type of cblk */
131  int8_t gpuid; /**< Store on which GPU the cblk is computed */
132  pastix_int_t fcolnum; /**< First column index (Global numbering) */
133  pastix_int_t lcolnum; /**< Last column index (Global numbering, inclusive) */
134  SolverBlok *fblokptr; /**< First block in column (diagonal) */
135  pastix_int_t stride; /**< Column block stride */
136  pastix_int_t lcolidx; /**< First column index (Local numbering), used for the rhs vectors */
137  pastix_int_t brownum; /**< First block in row facing the diagonal block in browtab, 0-based */
138  pastix_int_t brown2d; /**< First 2D-block in row facing the diagonal block in browtab, 0-based */
139  pastix_int_t sndeidx; /**< Global index of the original supernode the cblk belongs to */
140  pastix_int_t gcblknum; /**< Global column block index */
141  pastix_int_t bcscnum; /**< Index in the bcsctab if local cblk, -1 otherwise (FANIN | RECV) */
142  void *lcoeftab; /**< Coefficients access vector, lower part */
143  void *ucoeftab; /**< Coefficients access vector, upper part */
144  void *handler[2]; /**< Runtime data handler */
145  pastix_int_t selevtx; /**< Index to identify selected cblk for which intra-separator contributions are not compressed */
146  int ownerid; /**< Rank of the owner */
147  int threadid; /**< Rank of the accessing thread */
148 } SolverCblk;
149 
152 
155 
158 
159 /**
160  * @brief Solver column block structure.
161  *
162  * This structure stores all the numerical information about the factorization,
163  * as well as the structure of the problem. Only local information to each
164  * process is stored in this structure.
165  *
166  */
168  int restore; /**< Flag to indicate if it is require to restore data with
169  solverBackupRestore: 0: No need, 1:After solve,
170  2:After Factorization */
171  pastix_int_t baseval; /**< Base value for numberings */
172  pastix_int_t nodenbr; /**< Number of nodes before dof extension */
173  pastix_int_t coefnbr; /**< Number of coefficients (node after dof extension) */
174  pastix_int_t gcblknbr; /**< Global number of column blocks */
175  pastix_int_t cblknbr; /**< Local number of column blocks */
176  pastix_int_t faninnbr; /**< Local number of fanin cblk (included in cblknbr) */
177  pastix_int_t fanincnt; /**< Number of sends to realize */
178  pastix_int_t maxrecv; /**< Maximum blok size for a cblk_recv */
179  pastix_int_t recvnbr; /**< Local number of recv cblk (included in cblknbr) */
180  pastix_int_t recvcnt; /**< Number of receptions to realize */
181  pastix_int_t cblkmax1d; /**< Rank of the last cblk not beeing enabled for 2D computations */
182  pastix_int_t cblkmin2d; /**< Rank of the first cblk beeing enabled for 2D computations */
183  pastix_int_t cblkmaxblk; /**< Maximum number of blocks per cblk */
184  pastix_int_t cblkschur; /**< Index of the first local cblk in Schur */
185  pastix_int_t nb2dcblk; /**< Number of 2D cblks */
186  pastix_int_t nb2dblok; /**< Number of 2D blocks */
187  pastix_int_t bloknbr; /**< Number of blocks */
188  pastix_int_t brownbr; /**< Size of the browtab array */
189  SolverCblk * restrict cblktab; /**< Array of solver column blocks [+1] */
190  SolverBlok * restrict bloktab; /**< Array of solver blocks [+1] */
191  pastix_int_t * restrict browtab; /**< Array of blocks */
192 
193  pastix_int_t *gcbl2loc; /**< Array of local cblknum corresponding to gcblknum */
194 
195  pastix_lr_t lowrank; /**< Low-rank parameters */
196  pastix_factotype_t factotype; /**< General or symmetric factorization? */
197  double diagthreshold; /**< Diagonal threshold for pivoting */
198  volatile int32_t nbpivots; /**< Number of pivots during the factorization */
199 
200 #if defined(PASTIX_WITH_PARSEC)
201  parsec_sparse_matrix_desc_t *parsec_desc;
202 #endif
203 #if defined(PASTIX_WITH_STARPU)
204  starpu_sparse_matrix_desc_t *starpu_desc;
205  starpu_dense_matrix_desc_t *starpu_desc_rhs;
206 #endif
207 
208  pastix_int_t offdmax; /*+ Maximum size of the off-diagonal blocks for hetrf/sytrf temporary buffers +*/
209  pastix_int_t gemmmax; /*+ Maximum size of the GEMM update for 1d GEMM computations +*/
210  pastix_int_t blokmax; /*+ Maximum size of 2D blocks +*/
211  pastix_int_t colmax; /*+ Maximum column width in solvmtx +*/
212 
213  int clustnum; /*+ current processor number +*/
214  int clustnbr; /*+ number of processors +*/
215  pastix_int_t procnbr; /*+ Number of physical processor used +*/
216  pastix_int_t thrdnbr; /*+ Number of local computation threads +*/
217  pastix_int_t bublnbr; /*+ Number of local computation threads +*/
218  /* BubbleTree * restrict btree; /\*+ Bubbles tree +*\/ */
219 
220  Task * restrict tasktab; /*+ Task access vector +*/
221  pastix_int_t tasknbr; /*+ Number of Tasks +*/
222  pastix_int_t ** ttsktab; /*+ Task access vector by thread +*/
223  pastix_int_t * ttsknbr; /*+ Number of tasks by thread +*/
224  pastix_queue_t ** computeQueue; /*+ Queue of task to compute by thread +*/
225 
226  pastix_int_t *selevtx; /*+ Array to identify which cblk are pre-selected +*/
227 
228  MPI_Request *reqtab; /**< Array of requests for MPI asynchronous messages */
229  pastix_int_t *reqidx; /**< Array of local cblknum index corresponding to reqtab */
230  pastix_int_t reqnbr; /**< Length of the reqtab/reqidx arrays */
231  volatile int32_t reqnum; /**< Current amount of active requests (packed) */
232  pastix_atomic_lock_t reqlock; /**< Lock to access the request arrays */
233  void *rcoeftab; /**< Reception buffer for the communication */
234 
235  size_t *com_vector; /**< Matrix of communications between nodes. */
236 
237  PASTIX_Comm solv_comm; /*+ Copy of the pastix_data->inter_node_comm +*/
238 };
239 
240 /**
241  * @brief Compute the number of columns in a column block.
242  * @param[in] cblk
243  * The pointer to the column block.
244  * @return The number of columns in the cblk.
245  */
246 static inline pastix_int_t
247 cblk_colnbr( const SolverCblk *cblk )
248 {
249  return cblk->lcolnum - cblk->fcolnum + 1;
250 }
251 
252 /**
253  * @brief Get the pointer to the data associated to the lower part of the cblk.
254  * @param[in] cblk
255  * The pointer to the column block.
256  * @return lcoeftab if the cblk is not compressed.
257  * fblokptr->LRblock[0] if the cblk is compressed.
258  */
259 static inline void *
260 cblk_getdataL( const SolverCblk *cblk ) {
261  return (cblk->cblktype & CBLK_COMPRESSED) ? cblk->fblokptr->LRblock[0] : cblk->lcoeftab;
262 }
263 
264 /**
265  * @brief Get the pointer to the data associated to the upper part of the cblk.
266  * @param[in] cblk
267  * The pointer to the column block.
268  * @return lcoeftab if the cblk is not compressed.
269  * fblokptr->LRblock[1] if the cblk is compressed.
270  */
271 static inline void *
272 cblk_getdataU( const SolverCblk *cblk ) {
273  return (cblk->cblktype & CBLK_COMPRESSED) ? cblk->fblokptr->LRblock[1] : cblk->ucoeftab;
274 }
275 
276 /**
277  * @brief Get the pointer to the data associated to the side part of the cblk.
278  * @param[in] cblk
279  * The pointer to the column block.
280  * @param[in] side
281  * side of the data needed.
282  * PastixLCoef for the lower part PastixUCoef for the upper part.
283  * @return [lu]coeftab if the cblk is not compressed.
284  * fblokptr->LRblock[side] if the cblk is compressed.
285  */
286 static inline void *
288  return cblk->cblktype & CBLK_COMPRESSED
289  ? cblk->fblokptr->LRblock[side]
290  : side == PastixUCoef ? cblk->ucoeftab
291  : cblk->lcoeftab;
292 }
293 
294 /**
295  * @brief Compute the number of blocks in a column block.
296  * @param[in] cblk
297  * The pointer to the column block.
298  * @return The number of blocks in the cblk including the diagonal block.
299  */
300 static inline pastix_int_t
301 cblk_bloknbr( const SolverCblk *cblk )
302 {
303  return (cblk+1)->fblokptr - cblk->fblokptr + 1;
304 }
305 
306 /**
307  * @brief Compute the number of rows of a block.
308  * @param[in] blok
309  * The pointer to the block.
310  * @return The number of rows in the block.
311  */
312 static inline pastix_int_t
313 blok_rownbr( const SolverBlok *blok )
314 {
315  return blok->lrownum - blok->frownum + 1;
316 }
317 
318 /**
319  * @brief Compute the number of rows of a contiguous block in front of the same cblk.
320  * @param[in] blok
321  * The pointer to the block.
322  * @return The number of rows in the block.
323  */
324 static inline pastix_int_t
326 {
327  pastix_int_t rownbr = blok_rownbr( blok );
328 
329  while( (blok[0].lcblknm == blok[1].lcblknm) &&
330  (blok[0].fcblknm == blok[1].fcblknm) )
331  {
332  blok++;
333  rownbr += blok_rownbr( blok );
334  }
335  return rownbr;
336 }
337 
338 /**
339  * @brief Return if a block is preselected as either part of the projection, or
340  * as a sub-diagonal block.
341  * @param[in] cblk
342  * The pointer to the cblk to which belong the block.
343  * @param[in] blok
344  * The pointer to the block.
345  * @param[in] fcbk
346  * The pointer to the facing cblk of the block.
347  * @return True is the block is preselected, false otherwise.
348  */
349 static inline int
351  const SolverBlok *blok,
352  const SolverCblk *fcbk )
353 {
354  int is_preselected = ( fcbk->selevtx );
355  int is_firstoffd = ( blok == (cblk->fblokptr + 1) );
356 
357  return ( fcbk->sndeidx == cblk->sndeidx ) && ( is_preselected | is_firstoffd );
358 }
359 
360 /**
361  * @brief Compute the number of rows of a column block.
362  * @param[in] cblk
363  * The pointer to the column block.
364  * @return The number of rows in the column block.
365  */
366 static inline pastix_int_t
367 cblk_rownbr( const SolverCblk *cblk )
368 {
369  pastix_int_t rownbr = 0;
370  SolverBlok * blok;
371  for (blok = cblk->fblokptr; blok < cblk[1].fblokptr; blok++) {
372  rownbr += blok_rownbr(blok);
373  }
374  return rownbr;
375 }
376 
377 /**
378  * @brief Task stealing method.
379  *
380  * @param[inout] solvmtx
381  * The pointer to the solverMatrix.
382  * @param[in] rank
383  * Rank of the computeQueue.
384  * @param[inout] dest
385  * Rank of the stolen queue.
386  * @param[in] nbthreads
387  * Total amount of threads in the node.
388  * @return The concerned cblk if it exists. -1 otherwhise.
389  */
390 static inline pastix_int_t
391 stealQueue( SolverMatrix *solvmtx,
392  int rank,
393  int *dest,
394  int nbthreads )
395 {
396  int rk = *dest;
397  pastix_queue_t *stoleQueue;
398  pastix_int_t cblknum = -1;
399  while( rk != rank )
400  {
401  assert( solvmtx->computeQueue[ rk ] );
402  stoleQueue = solvmtx->computeQueue[ rk ];
403  if( (cblknum = pqueuePop(stoleQueue)) != -1 ){
404  *dest = rk;
405  return cblknum;
406  }
407  rk = (rk + 1)%nbthreads;
408  }
409  return cblknum;
410 }
411 
412 /**
413  * @brief Check if a block is included inside another one.
414  *
415  * Indicate if a blok is included inside an other block.
416  * i.e. indicate if the row range of the first block is included in the
417  * one of the second.
418  *
419  * @param[in] blok The block that is tested for inclusion.
420  * @param[in] fblok The block that is suppose to include the first one.
421  *
422  * @retval true if the first block is included in the second one.
423  * @retval false if the first block is not included in the second one.
424  */
425 static inline int
427  const SolverBlok *fblok )
428 {
429 # if defined(NAPA_SOPALIN)
430  return (((blok->frownum >= fblok->frownum) &&
431  (blok->lrownum <= fblok->lrownum)) ||
432  ((blok->frownum <= fblok->frownum) &&
433  (blok->lrownum >= fblok->lrownum)) ||
434  ((blok->frownum <= fblok->frownum) &&
435  (blok->lrownum >= fblok->frownum)) ||
436  ((blok->frownum <= fblok->lrownum) &&
437  (blok->lrownum >= fblok->lrownum)));
438 # else
439  return ((blok->frownum >= fblok->frownum) &&
440  (blok->lrownum <= fblok->lrownum));
441 # endif /* defined(NAPA_SOPALIN) */
442 }
443 
444 void solverInit( SolverMatrix *solvmtx );
445 void solverExit( SolverMatrix *solvmtx );
446 
447 int solverMatrixGen( SolverMatrix *solvmtx,
448  const symbol_matrix_t *symbmtx,
449  const pastix_order_t *ordeptr,
450  const SimuCtrl *simuctl,
451  const BlendCtrl *ctrl,
452  PASTIX_Comm comm,
453  isched_t *isched );
454 
455 int solverMatrixGenSeq( SolverMatrix *solvmtx,
456  const symbol_matrix_t *symbmtx,
457  const pastix_order_t *ordeptr,
458  const SimuCtrl *simuctl,
459  const BlendCtrl *ctrl,
460  PASTIX_Comm comm,
461  isched_t *isched,
462  pastix_int_t is_dbg );
463 
464 int solverLoad( SolverMatrix *solvptr,
465  FILE *stream );
466 int solverSave( const SolverMatrix *solvptr,
467  FILE *stream );
468 
469 void solverRealloc( SolverMatrix *solvptr);
470 SolverMatrix *solverCopy ( const SolverMatrix *solvptr,
471  int flttype );
472 
473 int solverCheck ( const SolverMatrix *solvmtx );
474 int solverDraw ( const SolverMatrix *solvptr,
475  FILE *stream,
476  int verbose,
477  const char *directory );
478 void solverPrintStats( const SolverMatrix *solvptr );
479 
480 void solverRequestInit( SolverMatrix *solvmtx );
481 void solverRequestExit( SolverMatrix *solvmtx );
482 
484  SolverMatrix *solvmtx,
485  int flttype );
486 void solverRecvExit( SolverMatrix *solvmtx );
487 
488 /*
489  * Solver backup
490  */
491 struct SolverBackup_s;
492 typedef struct SolverBackup_s SolverBackup_t;
493 
494 SolverBackup_t *solverBackupInit ( const SolverMatrix *solvmtx );
495 int solverBackupRestore( SolverMatrix *solvmtx, const SolverBackup_t *b );
496 void solverBackupExit ( SolverBackup_t *b );
497 
498 #endif /* _solver_h_ */
solver_cblk_s::ownerid
int ownerid
Definition: solver.h:146
solver_cblk_s::threadid
int threadid
Definition: solver.h:147
solverSave
int solverSave(const SolverMatrix *solvptr, FILE *stream)
Save a solver matrix structure into a file.
Definition: solver_io.c:261
solver_blok_s::frownum
pastix_int_t frownum
Definition: solver.h:112
cblk_rownbr
static pastix_int_t cblk_rownbr(const SolverCblk *cblk)
Compute the number of rows of a column block.
Definition: solver.h:367
solverBackupExit
void solverBackupExit(SolverBackup_t *b)
Free the solver backup data structure.
Definition: solver_backup.c:207
solver_matrix_s::browtab
pastix_int_t *restrict browtab
Definition: solver.h:191
symbol_matrix_s
Symbol matrix structure.
Definition: symbol.h:75
solver_matrix_s::cblkmaxblk
pastix_int_t cblkmaxblk
Definition: solver.h:183
blok_rownbr
static pastix_int_t blok_rownbr(const SolverBlok *blok)
Compute the number of rows of a block.
Definition: solver.h:313
cblk_colnbr
static pastix_int_t cblk_colnbr(const SolverCblk *cblk)
Compute the number of columns in a column block.
Definition: solver.h:247
SolverCblk
struct solver_cblk_s SolverCblk
Solver column block structure.
pastix_lr_s
Structure to define the type of function to use for the low-rank kernels and their parameters.
Definition: pastix_lowrank.h:147
solver_cblk_s::fblokptr
SolverBlok * fblokptr
Definition: solver.h:134
solverMatrixGenSeq
int solverMatrixGenSeq(SolverMatrix *solvmtx, const symbol_matrix_t *symbmtx, const pastix_order_t *ordeptr, const SimuCtrl *simuctl, const BlendCtrl *ctrl, PASTIX_Comm comm, isched_t *isched, pastix_int_t is_dbg)
Initialize the solver matrix structure in sequential.
Definition: solver_matrix_gen.c:424
pastix_queue_s
Queue structure.
Definition: queue.h:38
solver_blok_s::gpuid
int8_t gpuid
Definition: solver.h:116
pastix_lowrank.h
solver_blok_s::lcblknm
pastix_int_t lcblknm
Definition: solver.h:109
solverRequestExit
void solverRequestExit(SolverMatrix *solvmtx)
Free the arrays related to the requests.
Definition: solver.c:459
solver_cblk_s::stride
pastix_int_t stride
Definition: solver.h:135
solver_matrix_s::gcbl2loc
pastix_int_t * gcbl2loc
Definition: solver.h:193
solver_matrix_s::gcblknbr
pastix_int_t gcblknbr
Definition: solver.h:174
solver_matrix_s::recvcnt
pastix_int_t recvcnt
Definition: solver.h:180
solver_cblk_s
Solver column block structure.
Definition: solver.h:127
solverInit
void solverInit(SolverMatrix *solvmtx)
Initialize the solver structure.
Definition: solver.c:118
solver_blok_s::gbloknm
pastix_int_t gbloknm
Definition: solver.h:111
pastix_coefside_t
enum pastix_coefside_e pastix_coefside_t
Data blocks used in the kernel.
solverRecvExit
void solverRecvExit(SolverMatrix *solvmtx)
Free the array linked to pending reception.
Definition: solver.c:540
pastix_order_s
Order structure.
Definition: order.h:45
solver_cblk_recv_s
Solver recv column block structure.
Definition: solver.h:79
solver_blok_s::handler
void * handler[2]
Definition: solver.h:108
solver_matrix_s::reqtab
MPI_Request * reqtab
Definition: solver.h:228
solver_cblk_s::bcscnum
pastix_int_t bcscnum
Definition: solver.h:141
solver_blok_s
Solver block structure.
Definition: solver.h:107
solver_matrix_s::lowrank
pastix_lr_t lowrank
Definition: solver.h:195
solver_matrix_s::bloknbr
pastix_int_t bloknbr
Definition: solver.h:187
solver_matrix_s::nb2dcblk
pastix_int_t nb2dcblk
Definition: solver.h:185
solver_matrix_s::nodenbr
pastix_int_t nodenbr
Definition: solver.h:172
task_s::bloknum
pastix_int_t bloknum
Definition: solver.h:94
pastix_lrblock_s
The block low-rank structure to hold a matrix in low-rank form.
Definition: pastix_lowrank.h:112
solver_matrix_s::cblktab
SolverCblk *restrict cblktab
Definition: solver.h:189
task_s::ctrbcnt
pastix_int_t volatile ctrbcnt
Definition: solver.h:95
task_s::cblknum
pastix_int_t cblknum
Definition: solver.h:93
solver_cblk_s::ucoeftab
void * ucoeftab
Definition: solver.h:143
blendctrl_s
The type and structure definitions.
Definition: blendctrl.h:28
solver_matrix_s::rcoeftab
void * rcoeftab
Definition: solver.h:233
solver_cblk_recv_t
struct solver_cblk_recv_s solver_cblk_recv_t
Solver recv column block structure.
solver_cblk_recv_s::lcolnum
pastix_int_t lcolnum
Definition: solver.h:83
solver_blok_recv_s::frownum
pastix_int_t frownum
Definition: solver.h:70
solver_cblk_s::lcolnum
pastix_int_t lcolnum
Definition: solver.h:133
solverMatrixGen
int solverMatrixGen(SolverMatrix *solvmtx, const symbol_matrix_t *symbmtx, const pastix_order_t *ordeptr, const SimuCtrl *simuctl, const BlendCtrl *ctrl, PASTIX_Comm comm, isched_t *isched)
Initialize the solver matrix structure.
Definition: solver_matrix_gen.c:81
is_block_inside_fblock
static int is_block_inside_fblock(const SolverBlok *blok, const SolverBlok *fblok)
Check if a block is included inside another one.
Definition: solver.h:426
solver_matrix_s::reqnbr
pastix_int_t reqnbr
Definition: solver.h:230
solver_matrix_s::nbpivots
volatile int32_t nbpivots
Definition: solver.h:198
solver_blok_s::iluklvl
int iluklvl
Definition: solver.h:118
simuctrl_s
Control structure for the simulation.
Definition: simu.h:116
starpu_dense_matrix_desc_s
StarPU descriptor for the vectors linked to a given sparse matrix.
Definition: pastix_starpu.h:104
task_s
The task structure for the numerical factorization.
Definition: solver.h:90
solverPrintStats
void solverPrintStats(const SolverMatrix *solvptr)
Print statistical information about the solver matrix structure.
Definition: solver.c:196
solver_blok_s::browind
pastix_int_t browind
Definition: solver.h:115
blok_rownbr_ext
static pastix_int_t blok_rownbr_ext(const SolverBlok *blok)
Compute the number of rows of a contiguous block in front of the same cblk.
Definition: solver.h:325
solver_cblk_s::sndeidx
pastix_int_t sndeidx
Definition: solver.h:139
task_s::taskid
pastix_int_t taskid
Definition: solver.h:91
cblk_bloknbr
static pastix_int_t cblk_bloknbr(const SolverCblk *cblk)
Compute the number of blocks in a column block.
Definition: solver.h:301
solver_cblk_s::lcoeftab
void * lcoeftab
Definition: solver.h:142
solverDraw
int solverDraw(const SolverMatrix *solvptr, FILE *stream, int verbose, const char *directory)
Writes a PostScript picture of the low-rank solver matrix.
Definition: solver_draw.c:54
blok_is_preselected
static int blok_is_preselected(const SolverCblk *cblk, const SolverBlok *blok, const SolverCblk *fcbk)
Return if a block is preselected as either part of the projection, or as a sub-diagonal block.
Definition: solver.h:350
pastix_factotype_t
enum pastix_factotype_e pastix_factotype_t
Factorization algorithms available for IPARM_FACTORIZATION parameter.
solverRecvInit
void solverRecvInit(pastix_coefside_t side, SolverMatrix *solvmtx, int flttype)
Allocate the reception buffer, and initiate the first persistant reception.
Definition: solver.c:500
solver_cblk_recv_s::bloktab
solver_blok_recv_t bloktab[1]
Definition: solver.h:84
solverRequestInit
void solverRequestInit(SolverMatrix *solvmtx)
Instanciate the arrays for the requests according to the scheduler.
Definition: solver.c:423
PastixUCoef
@ PastixUCoef
Definition: api.h:457
solver_cblk_s::brown2d
pastix_int_t brown2d
Definition: solver.h:138
solver_matrix_s::nb2dblok
pastix_int_t nb2dblok
Definition: solver.h:186
solver_matrix_s::coefnbr
pastix_int_t coefnbr
Definition: solver.h:173
solver_matrix_s::maxrecv
pastix_int_t maxrecv
Definition: solver.h:178
solver_cblk_s::brownum
pastix_int_t brownum
Definition: solver.h:137
solverCheck
int solverCheck(const SolverMatrix *solvmtx)
Checks the consistency of the given solver matrix structure.
Definition: solver_check.c:54
solverExit
void solverExit(SolverMatrix *solvmtx)
Free the content of the solver matrix structure.
Definition: solver.c:143
solver_matrix_s::restore
int restore
Definition: solver.h:168
solver_cblk_recv_s::fcolnum
pastix_int_t fcolnum
Definition: solver.h:82
solver_blok_recv_s
Solver recv block structure.
Definition: solver.h:69
solver_cblk_s::ctrbcnt
volatile int32_t ctrbcnt
Definition: solver.h:129
solver_cblk_s::selevtx
pastix_int_t selevtx
Definition: solver.h:145
solver_matrix_s::bloktab
SolverBlok *restrict bloktab
Definition: solver.h:190
solver_matrix_s::brownbr
pastix_int_t brownbr
Definition: solver.h:188
solver_cblk_s::gcblknum
pastix_int_t gcblknum
Definition: solver.h:140
solver_blok_s::coefind
pastix_int_t coefind
Definition: solver.h:114
solver_blok_s::lrownum
pastix_int_t lrownum
Definition: solver.h:113
parsec_sparse_matrix_desc_s
Definition: pastix_parsec.h:33
cblk_getdata
static void * cblk_getdata(const SolverCblk *cblk, pastix_coefside_t side)
Get the pointer to the data associated to the side part of the cblk.
Definition: solver.h:287
solver_matrix_s::baseval
pastix_int_t baseval
Definition: solver.h:171
solver_matrix_s::cblkmin2d
pastix_int_t cblkmin2d
Definition: solver.h:182
Task
struct task_s Task
The task structure for the numerical factorization.
solver_matrix_s::reqnum
volatile int32_t reqnum
Definition: solver.h:231
solver_matrix_s::cblknbr
pastix_int_t cblknbr
Definition: solver.h:175
solver_cblk_s::cblktype
int8_t cblktype
Definition: solver.h:130
solver_cblk_s::lcolidx
pastix_int_t lcolidx
Definition: solver.h:136
task_s::prionum
pastix_int_t prionum
Definition: solver.h:92
solver_blok_s::inlast
int8_t inlast
Definition: solver.h:117
solver_cblk_s::lock
pastix_atomic_lock_t lock
Definition: solver.h:128
solver_cblk_s::fcolnum
pastix_int_t fcolnum
Definition: solver.h:132
solver_matrix_s::recvnbr
pastix_int_t recvnbr
Definition: solver.h:179
solver_matrix_s::cblkschur
pastix_int_t cblkschur
Definition: solver.h:184
solver_matrix_s::fanincnt
pastix_int_t fanincnt
Definition: solver.h:177
solverCopy
SolverMatrix * solverCopy(const SolverMatrix *solvptr, int flttype)
Generate a copy of a solver matrix structure.
Definition: solver_copy.c:171
solver_matrix_s
Solver column block structure.
Definition: solver.h:167
solver_blok_s::LRblock
pastix_lrblock_t * LRblock[2]
Definition: solver.h:121
solver_matrix_s::faninnbr
pastix_int_t faninnbr
Definition: solver.h:176
solverLoad
int solverLoad(SolverMatrix *solvptr, FILE *stream)
Load a solver matrix structure from a file.
Definition: solver_io.c:45
solver_matrix_s::com_vector
size_t * com_vector
Definition: solver.h:235
solverRealloc
void solverRealloc(SolverMatrix *solvptr)
Realloc in a contiguous way a given solver structure.
Definition: solver_copy.c:205
solver_matrix_s::diagthreshold
double diagthreshold
Definition: solver.h:197
solver_cblk_s::gpuid
int8_t gpuid
Definition: solver.h:131
solverBackupInit
SolverBackup_t * solverBackupInit(const SolverMatrix *solvmtx)
Initialize the backup structure.
Definition: solver_backup.c:59
pqueuePop
static pastix_int_t pqueuePop(pastix_queue_t *q)
Pop the head of the queue whithout returning the keys.
Definition: queue.h:75
solverBackupRestore
int solverBackupRestore(SolverMatrix *solvmtx, const SolverBackup_t *b)
Restore initial values.
Definition: solver_backup.c:140
solver_blok_recv_s::lrownum
pastix_int_t lrownum
Definition: solver.h:71
solver_matrix_s::factotype
pastix_factotype_t factotype
Definition: solver.h:196
solver_matrix_s::reqidx
pastix_int_t * reqidx
Definition: solver.h:229
solver_cblk_s::handler
void * handler[2]
Definition: solver.h:144
cblk_getdataU
static void * cblk_getdataU(const SolverCblk *cblk)
Get the pointer to the data associated to the upper part of the cblk.
Definition: solver.h:272
solver_matrix_s::reqlock
pastix_atomic_lock_t reqlock
Definition: solver.h:232
starpu_sparse_matrix_desc_s
StarPU descriptor stucture for the sparse matrix.
Definition: pastix_starpu.h:92
solver_matrix_s::cblkmax1d
pastix_int_t cblkmax1d
Definition: solver.h:181
cblk_getdataL
static void * cblk_getdataL(const SolverCblk *cblk)
Get the pointer to the data associated to the lower part of the cblk.
Definition: solver.h:260
solver_blok_s::fcblknm
pastix_int_t fcblknm
Definition: solver.h:110
stealQueue
static pastix_int_t stealQueue(SolverMatrix *solvmtx, int rank, int *dest, int nbthreads)
Task stealing method.
Definition: solver.h:391
SolverBlok
struct solver_blok_s SolverBlok
Solver block structure.
solver_blok_recv_t
struct solver_blok_recv_s solver_blok_recv_t
Solver recv block structure.