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