PaStiX Handbook  6.3.2
solver_matrix_gen_utils.c
Go to the documentation of this file.
1 /**
2  *
3  * @file solver_matrix_gen_utils.c
4  *
5  * PaStiX solver structure generation functions to factorize
6  * solver_matric_gen.c .
7  *
8  * @copyright 1998-2023 Bordeaux INP, CNRS (LaBRI UMR 5800), Inria,
9  * Univ. Bordeaux. All rights reserved.
10  *
11  * @version 6.3.2
12  * @author Tony Delarue
13  * @author Pascal Henon
14  * @author Pierre Ramet
15  * @author Xavier Lacoste
16  * @author Mathieu Faverge
17  * @author Alycia Lisito
18  * @author Nolan Bredel
19  * @date 2023-12-01
20  *
21  * @addtogroup blend_dev_solver
22  * @{
23  *
24  **/
25 #include "common.h"
26 #include "symbol/symbol.h"
27 #include "blend/solver.h"
28 #include "elimintree.h"
29 #include "cost.h"
30 #include "cand.h"
31 #include "pastix/order.h"
32 #include "extendVector.h"
33 #include "simu.h"
35 
36 /**
37  *******************************************************************************
38  *
39  * @brief Fill the local numbering arrays to compress the symbol information
40  * into solver.
41  *
42  *******************************************************************************
43  *
44  * @param[in] symbmtx
45  * The pointer to the symbol matrix structure.
46  *
47  * @param[in] simuctrl
48  * The pointer to the simuctrl structure.
49  *
50  * @param[inout] solvmtx
51  * Pointer to the solver matrix.
52  *
53  * @param[inout] cblklocalnum
54  * Local cblk infos.
55  *
56  * @param[inout] bloklocalnum
57  * Local blok infos.
58  *
59  * @param[inout] tasklocalnum
60  * Local tasks infos.
61  *
62  * @param[inout] ftgttab
63  * Array of fan-in to store the lists of recv/fanin cblk per local cblk.
64  *
65  *******************************************************************************/
66 void
68  const SimuCtrl *simuctrl,
69  SolverMatrix *solvmtx,
70  pastix_int_t *cblklocalnum,
71  pastix_int_t *bloklocalnum,
72  pastix_int_t *tasklocalnum,
73  solver_cblk_recv_t **ftgttab )
74 {
75  pastix_int_t *localindex;
76  symbol_cblk_t *symbcblk;
77  symbol_blok_t *symbblok;
78  pastix_int_t cblknum, brownum, brownbr;
79  pastix_int_t faninnbr, recvnbr;
80  pastix_int_t i, j, k, c, fc;
81  pastix_int_t flaglocal;
82  pastix_int_t clustnum = solvmtx->clustnum;
83 
84  /* Initialize the set of cluster candidates for each cblk */
85  MALLOC_INTERN( localindex, solvmtx->clustnbr, pastix_int_t );
86 
87  /*
88  * Compute local number of tasks on each cluster
89  */
90  memset( localindex, 0, solvmtx->clustnbr * sizeof(pastix_int_t) );
91  for ( i = 0; i < simuctrl->tasknbr; i++ ) {
92  c = simuctrl->bloktab[ simuctrl->tasktab[i].bloknum ].ownerclust;
93 
94  tasklocalnum[i] = localindex[c];
95  localindex[c]++;
96  }
97  solvmtx->tasknbr = localindex[clustnum];
98 
99  /*
100  * Compute the local numbering of the fan-in and recv cblks on each cluster
101  */
102  /* Reset the array to compute local informations */
103  memset( localindex, 0, solvmtx->clustnbr * sizeof( pastix_int_t ) );
104 
105  cblknum = 0;
106  brownum = 0;
107  recvnbr = 0;
108  faninnbr = 0;
109  symbcblk = symbmtx->cblktab;
110  for ( i = 0; i < symbmtx->cblknbr; i++, symbcblk++ ) {
111  brownbr = symbcblk[1].brownum - symbcblk[0].brownum;
112 
113  /*
114  * The cblk is considered local if data are local, or if we store a
115  * compressed copy for fanin
116  */
117  flaglocal = ( simuctrl->cblktab[i].owned ) || ( ftgttab[i] != NULL );
118  if ( !flaglocal ) {
119  cblklocalnum[i] = -i - 1;
120 #if !defined(NDEBUG)
121  for ( j=symbcblk[0].bloknum; j<symbcblk[1].bloknum; j++ ) {
122  bloklocalnum[j] = -1;
123  assert( simuctrl->bloktab[j].ownerclust != clustnum );
124  }
125 #endif
126  continue;
127  }
128 
129  /*
130  * The cblk is local.
131  */
132  if ( simuctrl->cblktab[i].owned ) {
133  solver_cblk_recv_t *ftgtcblk;
134 
135  /*
136  * The cblk is local. We may receive remote information, let's:
137  * - compute the size of the compressed browtab
138  * - compute the number of remote fanin to be received for the update
139  * - compute the set of remote nodes sending a fanin
140  *
141  * To do that, we work on the incoming edges.
142  */
143  for ( j = symbcblk[0].brownum; j < symbcblk[1].brownum; j++ ) {
144  k = symbmtx->browtab[j];
145  symbblok = symbmtx->bloktab + k;
146  c = simuctrl->bloktab[k].ownerclust;
147 
148  assert( i == symbblok->fcblknm );
149 
150  /* This is a remote contribution we add it to the ftgt and update the counters */
151  if ( c != clustnum ) {
152  solver_recv_update_recv( ftgttab + i,
153  symbmtx,
154  symbmtx->cblktab + symbblok->lcblknm,
155  symbblok, symbcblk, c );
156  brownbr--;
157  }
158  assert( brownbr >= 0 );
159  }
160 
161  /*
162  * Now that all remote contributions have been computed and summarized in ftgttab[i],
163  * let's compute the local information for the indices
164  */
165  ftgtcblk = ftgttab[i];
166  while( ftgtcblk != NULL ) {
167  assert( (ftgtcblk->ownerid != -1) &&
168  (ftgtcblk->ownerid != clustnum) );
169 
170  /* Book some space for the reception blocks */
171  localindex[clustnum] += solver_recv_get_bloknbr( ftgtcblk, symbcblk,
172  symbmtx->bloktab + symbcblk->bloknum );
173 
174  brownbr++; /* One more blok will be in the browtab */
175  cblknum++; /* Add one cblk */
176  recvnbr++; /* Add one reception count */
177 
178  ftgtcblk = ftgtcblk->next;
179  }
180 
181  /*
182  * Now, we need to get through the outgoing dependencies to generate
183  * the fanin informations if it needs to be added, and to compute
184  * local block indices.
185  */
186  symbblok = symbmtx->bloktab + symbcblk->bloknum;
187  for ( j=symbcblk[0].bloknum; j<symbcblk[1].bloknum; j++, symbblok++ ) {
188  symbol_cblk_t *symbfcbk;
189  pastix_int_t fcblknum, fbloknum;
190 
191  bloklocalnum[j] = localindex[clustnum];
192  localindex[clustnum]++;
193 
194  assert( simuctrl->bloktab[j].ownerclust == clustnum );
195 
196  fcblknum = symbblok->fcblknm;
197  symbfcbk = symbmtx->cblktab + fcblknum;
198  fbloknum = symbfcbk->bloknum;
199  fc = simuctrl->bloktab[fbloknum].ownerclust;
200 
201  /*
202  * If the facing cblk isn't local, we need to have a local copy
203  * of it to store the fan-in
204  */
205  if ( fc != clustnum ) {
206  solver_recv_update_fanin( ftgttab + fcblknum,
207  symbmtx, symbcblk, symbblok, symbfcbk, fc );
208  }
209  }
210  }
211  else {
212  /* If the cblk is not local, it is a fanin */
213 
214  /*
215  * First, let's look at incoming dependencies to reduce the brownbr
216  */
217  {
218  for ( j = symbcblk[0].brownum; j < symbcblk[1].brownum; j++ ) {
219  k = symbmtx->browtab[j];
220  c = simuctrl->bloktab[k].ownerclust;
221  if ( c != clustnum ) {
222  brownbr--;
223  }
224  }
225  }
226 
227  /*
228  * Second, let's update the localindex counter based on the number of local blocks
229  */
230  {
231  /* Check we have one and only one solver_cblk_recv associated to it */
232  solver_cblk_recv_t *ftgtcblk = ftgttab[i];
233  solver_blok_recv_t *ftgtblok;
234  assert( ftgtcblk != NULL );
235  assert( ftgtcblk->next == NULL );
236 
237  symbblok = symbmtx->bloktab + symbcblk->bloknum;
238  faninnbr++;
239  ftgtblok = ftgtcblk->bloktab;
240  for ( j=symbcblk[0].bloknum; j<symbcblk[1].bloknum; j++, symbblok++, ftgtblok++ )
241  {
242  assert( simuctrl->bloktab[j].ownerclust != clustnum );
243 
244  if ( (ftgtblok->frownum <= ftgtblok->lrownum) &&
245  (ftgtblok->frownum >= symbblok->frownum) &&
246  (ftgtblok->lrownum <= symbblok->lrownum) )
247  {
248  bloklocalnum[j] = localindex[clustnum];
249  localindex[clustnum]++;
250  }
251  else {
252  bloklocalnum[j] = -1;
253  }
254  }
255  }
256  }
257 
258  /* Store index of the current cblk */
259  cblklocalnum[i] = cblknum;
260  cblknum++;
261 
262  /* Update the brownum index */
263  brownum += brownbr;
264  assert( brownum <= symbcblk[1].brownum );
265  }
266 
267  solvmtx->cblknbr = cblknum;
268  solvmtx->bloknbr = localindex[clustnum];
269  solvmtx->brownbr = brownum;
270 
271  /* Reallocate recv_sources tab to diminish it's size */
272  solvmtx->recvnbr = recvnbr;
273  solvmtx->faninnbr = faninnbr;
274 
275  memFree_null( localindex );
276 }
277 
278 /**
279  *******************************************************************************
280  *
281  * @brief Reorder the browtab from the symbol structure in a distributed way.
282  * First stock the 1D blocks and then the 2D blocks.
283  *
284  *******************************************************************************
285  *
286  * @param[in] symbmtx
287  * The pointer to the symbol matrix structure.
288  *
289  * @param[in] symbcblk
290  * The pointer to the current symbol cblk.
291  *
292  * @param[inout] solvmtx
293  * Pointer to the solver matrix.
294  *
295  * @param[inout] solvcblk
296  * The pointer to the current solver cblk.
297  *
298  * @param[inout] browtmp
299  * Workspace array used to reorder the local brow information.
300  * Must be of size at least (symbcblk[1].brownum - symbcblk[0].brownum)
301  *
302  * @param[in] cblklocalnum
303  * Local cblk indices.
304  *
305  * @param[in] bloklocalnum
306  * Local blok indices.
307  *
308  * @param[in] brownum
309  * Current brownum.
310  *
311  *******************************************************************************
312  *
313  * @retval TODO
314  *
315  *******************************************************************************/
318  const symbol_cblk_t *symbcblk,
319  SolverMatrix *solvmtx,
320  SolverCblk *solvcblk,
321  pastix_int_t *browtmp,
322  const pastix_int_t *cblklocalnum,
323  const pastix_int_t *bloklocalnum,
324  pastix_int_t brownum )
325 {
326  pastix_int_t brownbr;
327  symbol_blok_t *symbblok;
328  SolverBlok *solvblok;
329  SolverCblk *browcblk;
330  pastix_int_t lcblknm, lbloknm;
331  pastix_int_t j2d, j1d, j, jmax;
332  pastix_int_t *b;
333 
334  brownbr = symbcblk[1].brownum - symbcblk[0].brownum;
335  solvcblk->brown2d = solvcblk->brownum + brownbr;
336 
337  /* Nothing to do here */
338  if ( !brownbr ) {
339  return 0;
340  }
341 
342  assert( brownbr <= symbmtx->browmax );
343  memcpy( browtmp,
344  symbmtx->browtab + symbcblk->brownum,
345  brownbr * sizeof(pastix_int_t) );
346 
347  /*
348  * j is the index in the local browtab (~ postition of b)
349  * j1d is the number of discovered 1d block in the browtab
350  * j2d if the index of the first 2d block in the original tab
351  * jmax is equal to brownbr
352  * brownbr is updated to store the real number of brow (minus fanin)
353  *
354  * b is a pointer to the temporary copy of the subsection of the browtab
355  * At the end of the first pass, if b[i] is negative, it has already been treated
356  * with if equal to:
357  * -1, the block was a 1D block
358  * -2, the block belonged to a remote cblk
359  * -3, the block belonged to a local fanin (should not happen)
360  * It is is positive, it's a 2D block that need to be pushed to the end of
361  * the browtab in the second pass.
362  */
363  b = browtmp;
364  j2d = -1;
365  jmax = brownbr;
366 
367  /* First pass to copy 1D updates */
368  for ( j=0, j1d=0; j < jmax; j++, b++ ) {
369  /* Get the contributing block in the symbol */
370  symbblok = symbmtx->bloktab + (*b);
371 
372  lcblknm = ( cblklocalnum == NULL ) ? symbblok->lcblknm : cblklocalnum[ symbblok->lcblknm ];
373 
374  /* If distant blok */
375  if ( lcblknm < 0 ) {
376  *b = -2;
377  brownbr--;
378  continue;
379  }
380 
381  /* Get the local cblk which owns the block */
382  browcblk = solvmtx->cblktab + lcblknm;
383 
384  /* Recv should never appear through cblklocalnum */
385  assert( !(browcblk->cblktype & CBLK_RECV) );
386 
387  /* Fanin should not contribute to local data */
388  if( browcblk->cblktype & CBLK_FANIN ) {
389  *b = -3;
390  brownbr--;
391  continue;
392  }
393 
394  /* Store the first non 1D index to not rediscover the begining, and skip 2d for now */
395  if ( browcblk->cblktype & CBLK_TASKS_2D ) {
396  j2d = ( j2d == -1 ) ? j : j2d;
397  continue;
398  }
399 
400  /* Find the SolvBlok corresponding to the SymbBlok */
401  lbloknm = ( bloklocalnum == NULL ) ? *b : bloklocalnum[ *b ];
402  solvblok = solvmtx->bloktab + lbloknm;
403 
404  assert( solvblok->lcblknm == lcblknm );
405 #if !defined(NDEBUG)
406  {
407  pastix_int_t frownum, lrownum;
408  symbol_blok_get_rownum( symbmtx, symbblok, &frownum, &lrownum );
409  assert( ( frownum == solvblok->frownum ) &&
410  ( lrownum == solvblok->lrownum ) );
411  }
412 #endif
413 
414  assert( brownum + j1d < solvmtx->brownbr );
415  solvmtx->browtab[brownum + j1d] = lbloknm;
416  solvblok->browind = brownum + j1d;
417  *b = -1;
418  j1d++;
419  }
420 
421  /* Store the index of the first 2D contribution in the array */
422  assert( j1d <= brownbr );
423  solvcblk->brown2d = solvcblk->brownum + j1d;
424 
425  /* Second pass to copy 2D updates to the end */
426  if ( j2d != -1 ) {
427  b = browtmp + j2d;
428  for ( j = j2d; j < jmax; j++, b++ ) {
429  symbblok = symbmtx->bloktab + ( *b );
430 
431  if ( *b < 0 ) {
432  continue;
433  }
434  lcblknm = ( cblklocalnum == NULL ) ? symbblok->lcblknm : cblklocalnum[ symbblok->lcblknm ];
435  assert( lcblknm >= 0 );
436 
437  /* Get the local cblk which owns the block */
438  browcblk = solvmtx->cblktab + lcblknm;
439  assert( (cblklocalnum == NULL) ||
440  (browcblk->ownerid == solvmtx->clustnum) );
441 
442  /* Find the SolvBlok corresponding to the SymbBlok */
443  lbloknm = ( bloklocalnum == NULL ) ? *b : bloklocalnum[ *b ];
444  solvblok = solvmtx->bloktab + lbloknm;
445 
446  assert( solvblok->lcblknm == lcblknm );
447  assert( ( symbblok->frownum == solvblok->frownum ) &&
448  ( symbblok->lrownum == solvblok->lrownum ) );
449 
450  assert( brownum + j1d < solvmtx->brownbr );
451  solvmtx->browtab[brownum + j1d] = lbloknm;
452  solvblok->browind = brownum + j1d;
453  j1d++;
454  }
455  }
456  assert( j1d == brownbr );
457 
458  return brownbr;
459 }
460 
461 /**
462  * @brief Structure to pass information to the muti-threaded ttsktab
463  * initialization function.
464  */
465 struct args_ttsktab
466 {
467  SolverMatrix *solvmtx; /**< Pointer to the solver matrix */
468  const SimuCtrl *simuctrl; /**< Pointer to simulation control structure */
469  const pastix_int_t *tasklocalnum; /**< Array of the local indices of the tasks */
470  pastix_int_t clustnum; /**< Index of the local cluster */
471 };
472 
473 /**
474  *******************************************************************************
475  *
476  * @brief Fill the ttsktab for it's own thread.
477  *
478  *******************************************************************************
479  *
480  * @param[in] ctx
481  * The context of the current thread
482  *
483  * @param[inout] args
484  * The pointer to the args_ttsktab structure that parameterize the
485  * function call.
486  *
487  *******************************************************************************/
488 void
489 solvMatGen_fill_ttsktab( isched_thread_t *ctx, void *args )
490 {
491  struct args_ttsktab *arg = (struct args_ttsktab*)args;
492  SolverMatrix *solvmtx = arg->solvmtx;
493  const SimuCtrl *simuctrl = arg->simuctrl;
494  const pastix_int_t *tasklocalnum = arg->tasklocalnum;
495  pastix_int_t clustnum = arg->clustnum;
496  int rank = ctx->rank;
497  SimuProc *simuproc = simuctrl->proctab
498  + ( simuctrl->clustab[clustnum].fprocnum + rank );
499  pastix_int_t i;
500  pastix_int_t priomin = PASTIX_INT_MAX;
501  pastix_int_t priomax = 0;
502  pastix_int_t ttsknbr = extendint_Size( simuproc->tasktab );
503  pastix_int_t j, jloc;
504 
505  solvmtx->ttsknbr[rank] = ttsknbr;
506  if(ttsknbr > 0) {
507  MALLOC_INTERN(solvmtx->ttsktab[rank], ttsknbr, pastix_int_t);
508  }
509  else {
510  solvmtx->ttsktab[rank] = NULL;
511  }
512 
513  for(i=0; i<ttsknbr; i++)
514  {
515  j = extendint_Read(simuproc->tasktab, i);
516  if( tasklocalnum != NULL ){
517  jloc = tasklocalnum[j];
518  }
519  else {
520  jloc = j;
521  }
522  /* Only local cblks should appear in the tasktab */
523  assert( !(solvmtx->cblktab[ solvmtx->tasktab[jloc].cblknum ].cblktype & (CBLK_FANIN|CBLK_RECV)) );
524  solvmtx->ttsktab[rank][i] = jloc;
525  solvmtx->cblktab[jloc].threadid = rank;
526 
527 #if defined(PASTIX_DYNSCHED)
528  solvmtx->tasktab[jloc].threadid = rank;
529 #endif
530  priomax = pastix_imax( solvmtx->tasktab[jloc].prionum, priomax );
531  priomin = pastix_imin( solvmtx->tasktab[jloc].prionum, priomin );
532  }
533 
534 #if defined(PASTIX_DYNSCHED)
535  solvmtx->btree->nodetab[rank].priomin = priomin;
536  solvmtx->btree->nodetab[rank].priomax = priomax;
537 #endif
538 }
539 
540 /**
541  *******************************************************************************
542  *
543  * @brief Fill in ttsktab for it's own thread. Only for debugging factorization.
544  *
545  *******************************************************************************
546  *
547  * @param[in] ctx
548  * the context of the current thread
549  *
550  * @param[inout] args
551  * The pointer to the args_ttsktab structure that parameterize the
552  * function call.
553  *
554  *******************************************************************************/
555 void
556 solvMatGen_fill_ttsktab_dbg( isched_thread_t *ctx, void *args )
557 {
558  struct args_ttsktab *arg = (struct args_ttsktab*)args;
559 
560  pastix_int_t i, j, size;
561  SolverMatrix *solvmtx = arg->solvmtx;
562  int rank = ctx->rank;
563  int nthread = ctx->global_ctx->world_size;
564  pastix_int_t tasknbr = solvmtx->tasknbr / nthread;
565  pastix_int_t priomin = PASTIX_INT_MAX;
566  pastix_int_t priomax = 0;
567 
568  size = (rank == nthread-1) ? (solvmtx->tasknbr - (nthread-1) * tasknbr) : tasknbr;
569  solvmtx->ttsknbr[rank] = size;
570 
571  if(size > 0) {
572  MALLOC_INTERN(solvmtx->ttsktab[rank], size, pastix_int_t);
573  }
574  else {
575  solvmtx->ttsktab[rank] = NULL;
576  }
577 
578  j = ((solvmtx->tasknbr - (nthread-1) * tasknbr) * rank);
579  for(i=0; i < size; i++)
580  {
581  solvmtx->ttsktab[rank][i] = j;
582 
583 #if defined(PASTIX_DYNSCHED)
584  solvmtx->tasktab[j].threadid = rank;
585 #endif
586  priomax = pastix_imax( solvmtx->tasktab[j].prionum, priomax );
587  priomin = pastix_imin( solvmtx->tasktab[j].prionum, priomin );
588  j++;
589  }
590 
591 #if defined(PASTIX_DYNSCHED)
592  solvmtx->btree->nodetab[rank].priomin = priomin;
593  solvmtx->btree->nodetab[rank].priomax = priomax;
594 #endif
595 }
596 
597 /**
598  *******************************************************************************
599  *
600  * @brief Fill the global tasktab array, as well as the thread ttsktab arrays.
601  *
602  *******************************************************************************
603  *
604  * @param[inout] solvmtx
605  * Pointer to the solver matrix.
606  *
607  * @param[in] isched
608  * The internal context to run multi-threaded functions.
609  *
610  * @param[in] simuctrl
611  * The pointer to the simulation control structure.
612  *
613  * @param[in] tasklocalnum
614  * Array of the local indices of the tasks.
615  *
616  * @param[in] cblklocalnum
617  * Array of the local indices of the cblk.
618  *
619  * @param[in] bloklocalnum
620  * Array of the local indices of the blocks.
621  *
622  * @param[in] clustnum
623  * Rank of the MPI instance.
624  *
625  * @param[in] is_dbg
626  * Enable/disable the ttsktab debug generation.
627  *
628  *******************************************************************************/
629 void
631  isched_t *isched,
632  const SimuCtrl *simuctrl,
633  const pastix_int_t *tasklocalnum,
634  const pastix_int_t *cblklocalnum,
635  const pastix_int_t *bloklocalnum,
636  pastix_int_t clustnum,
637  int is_dbg )
638 {
639  Task *solvtask;
640  SolverBlok *blok, *lblk;
641  SimuTask *simutask = simuctrl->tasktab;
642  SolverCblk *solvcblk = solvmtx->cblktab;
643  pastix_int_t tasknum = 0;
644  pastix_int_t tasknbr_1dp = 0;
645  pastix_int_t i;
646 
647  MALLOC_INTERN( solvmtx->tasktab, solvmtx->tasknbr+1, Task );
648  solvtask = solvmtx->tasktab;
649 
650  /* No local indices, this is a global solver */
651  if ( tasklocalnum == NULL )
652  {
653  for(i=0; i<simuctrl->tasknbr; i++, simutask++)
654  {
655  solvtask->taskid = COMP_1D;
656  solvtask->prionum = simutask->prionum;
657  solvtask->cblknum = simutask->cblknum;
658  solvtask->bloknum = simutask->bloknum;
659  solvtask->ctrbcnt = simutask->ctrbcnt;
660 
661  solvcblk = solvmtx->cblktab + solvtask->cblknum;
662  solvcblk->priority = solvtask->prionum;
663 
664  /* Schur, fanin and recv are counted only in static */
665  if ( solvcblk->cblktype & (CBLK_IN_SCHUR | CBLK_FANIN | CBLK_RECV ) ) {
666  tasknum++; solvtask++;
667  continue;
668  }
669 
670  /* Count the number of additional tasks in 1D+ mode */
671  if ( solvcblk->cblktype & CBLK_TASKS_2D )
672  {
673  blok = solvcblk[0].fblokptr + 1; /* first off-diagonal block */
674  lblk = solvcblk[1].fblokptr; /* the next diagonal block */
675 
676  /* if there are off-diagonal supernodes in the column */
677  for( ; blok < lblk; blok++ ) {
678  tasknbr_1dp++;
679  /* Skip blocks facing the same cblk */
680  while ( ( blok < lblk ) &&
681  ( blok[0].fcblknm == blok[1].fcblknm ) &&
682  ( blok[0].lcblknm == blok[1].lcblknm ) )
683  {
684  blok++;
685  }
686  }
687  }
688  tasknbr_1dp++;
689  tasknum++; solvtask++;
690  }
691  }
692  else
693  {
694  for(i=0; i<simuctrl->tasknbr; i++, simutask++)
695  {
696  if( simuctrl->bloktab[ simutask->bloknum ].ownerclust != clustnum )
697  {
698  continue;
699  }
700 
701  assert( tasknum == tasklocalnum[i] );
702 
703  solvtask->taskid = COMP_1D;
704  solvtask->prionum = simutask->prionum;
705  solvtask->cblknum = cblklocalnum[ simutask->cblknum ];
706  solvtask->bloknum = bloklocalnum[ simutask->bloknum ];
707  solvtask->ctrbcnt = simutask->ctrbcnt;
708 
709  solvcblk = solvmtx->cblktab + solvtask->cblknum;
710  solvcblk->priority = solvtask->prionum;
711 
712  /* Schur, fanin and recv are counted only in static */
713  if ( solvcblk->cblktype & (CBLK_IN_SCHUR | CBLK_FANIN | CBLK_RECV ) ) {
714  tasknum++; solvtask++;
715  continue;
716  }
717 
718  /* Count the number of additional tasks in 1D+ mode */
719  if ( solvcblk->cblktype & CBLK_TASKS_2D )
720  {
721  blok = solvcblk[0].fblokptr + 1; /* first off-diagonal block */
722  lblk = solvcblk[1].fblokptr; /* the next diagonal block */
723 
724  /* if there are off-diagonal supernodes in the column */
725  for( ; blok < lblk; blok++ ) {
726  tasknbr_1dp++;
727  /* Skip blocks facing the same cblk */
728  while ( ( blok < lblk ) &&
729  ( blok[0].fcblknm == blok[1].fcblknm ) &&
730  ( blok[0].lcblknm == blok[1].lcblknm ) )
731  {
732  blok++;
733  }
734  }
735  }
736 
737  tasknbr_1dp++;
738  tasknum++; solvtask++;
739  }
740  }
741  assert(tasknum == solvmtx->tasknbr);
742  solvmtx->tasknbr_1dp = tasknbr_1dp;
743 
744  /* One more task to avoid side effect */
745  solvtask->taskid = -1;
746  solvtask->prionum = -1;
747  solvtask->cblknum = solvmtx->cblknbr+1;
748  solvtask->bloknum = solvmtx->bloknbr+1;
749  solvtask->ctrbcnt = 0;
750 
751  /* Fill in the ttsktab arrays (one per thread) */
752  MALLOC_INTERN(solvmtx->ttsknbr, solvmtx->bublnbr, pastix_int_t );
753  MALLOC_INTERN(solvmtx->ttsktab, solvmtx->bublnbr, pastix_int_t* );
754 
755  if( is_dbg ) {
756  struct args_ttsktab args = { solvmtx, NULL, tasklocalnum, clustnum };
757  isched_parallel_call( isched, solvMatGen_fill_ttsktab_dbg, &args );
758  }
759  else {
760  struct args_ttsktab args = { solvmtx, simuctrl, tasklocalnum, clustnum };
761  isched_parallel_call( isched, solvMatGen_fill_ttsktab, &args );
762  }
763 }
764 
765 /**
766  *******************************************************************************
767  *
768  * @brief Compute the maximum area of the temporary buffers
769  * used during computation
770  *
771  * During this loop, we compute the maximum area that will be used as
772  * temporary buffers, and statistics:
773  * - diagmax: Only for hetrf/sytrf factorization, this the maximum size
774  * of a panel of MAXSIZEOFBLOCKS width in a diagonal block
775  * - gemmmax: For all, this is the maximum area used to compute the
776  * compacted gemm on a CPU.
777  *
778  * Rk: This loop is not merged within the main block loop, since strides have
779  * to be peviously computed.
780  *
781  *******************************************************************************
782  *
783  * @param[inout] solvmtx
784  * Pointer to the solver matrix.
785  *
786  *******************************************************************************/
787 void
789 {
790  SolverCblk *solvcblk = solvmtx->cblktab;
791  SolverBlok *solvblok = solvmtx->bloktab;
792  pastix_int_t gemmmax = 0;
793  pastix_int_t offdmax = 0;
794  pastix_int_t blokmax = 0;
795  pastix_int_t gemmarea, offdarea, cblk_m, acc_m, i;
796 
797  for(i=0; i<solvmtx->cblknbr; i++, solvcblk++)
798  {
799  SolverBlok *lblok = solvcblk[1].fblokptr;
800  pastix_int_t m = solvcblk->stride;
801  pastix_int_t n = cblk_colnbr( solvcblk );
802  pastix_int_t k = blok_rownbr( solvblok );
803 
804  m -= n;
805 
806  /*
807  * Compute the surface of the off-diagonal block in a panel for
808  * LDL^[th] factorizations
809  */
810  offdarea = m * n;
811  offdmax = pastix_imax( offdmax, offdarea );
812 
813  /*
814  * Compute the maximum area for 1d temporary workspace in GEMM
815  */
816  solvblok++;
817  cblk_m = -1;
818  acc_m = 0;
819  for( ; solvblok<lblok; solvblok++ ) {
820  k = blok_rownbr( solvblok );
821 
822  /*
823  * Temporary workspace for GEMM
824  * m+1 to store the diagonal in case of GEMDM
825  */
826  if ( !(solvcblk->cblktype & CBLK_LAYOUT_2D) ) {
827  gemmarea = (m+1) * k;
828  gemmmax = pastix_imax( gemmmax, gemmarea );
829  }
830 
831  /*
832  * Max size for off-diagonal blocks for 2-terms version of the
833  * 2D LDL
834  */
835  if ( solvcblk->cblktype & (CBLK_TASKS_2D | CBLK_COMPRESSED) ) {
836  if ( solvblok->fcblknm == cblk_m ) {
837  acc_m += k;
838  }
839  else {
840  cblk_m = solvblok->fcblknm;
841  acc_m = k;
842  }
843  /*
844  * acc_m+1 to store the diagonal in case of GEMDM
845  */
846  blokmax = pastix_imax( n * (acc_m+1), blokmax );
847  }
848  m -= k;
849  }
850  }
851 
852  solvmtx->offdmax = offdmax;
853  solvmtx->gemmmax = gemmmax;
854  solvmtx->blokmax = blokmax;
855 }
856 
857 /**
858  *******************************************************************************
859  *
860  * @brief Mark blocks if they belong to the last supernode, or if they are
861  * facing it for statistical purpose only.
862  *
863  * TODO : Should be improved by using the brow array in order to cover only the
864  * blocks in front of the last cblk
865  *
866  *******************************************************************************
867  *
868  * @param[inout] solvmtx
869  * Pointer to the solver matrix.
870  *
871  *******************************************************************************/
872 void
874 {
875 #if defined(PASTIX_SUPERNODE_STATS)
876  pastix_int_t i;
877  SolverBlok *solvblok = solvmtx->bloktab;
878 
879  for(i=0; i<solvmtx->bloknbr; i++, solvblok++ ) {
880  SolverCblk *fcblk = solvmtx->cblktab + solvblok->fcblknm;
881  SolverCblk *lcblk = solvmtx->cblktab + solvblok->lcblknm;
882  if ( fcblk->cblktype & CBLK_IN_LAST ) {
883  if ( lcblk->cblktype & CBLK_IN_LAST ) {
884  solvblok->inlast = 2;
885  }
886  else {
887  solvblok->inlast = 1;
888  }
889  }
890  }
891 #else
892  (void)solvmtx;
893 #endif
894 }
895 
896 /**
897  *******************************************************************************
898  *
899  * @brief Register a local cblk from a symbol_cblk_t structure !(Fanin|Recv)
900  *
901  *******************************************************************************
902  *
903  * @param[in] symbmtx
904  * The pointer to the symbol matrix.
905  *
906  * @param[in] candcblk
907  * The cand structure associated to the current cblk to get the type of
908  * the cblk.
909  *
910  * @param[in] cblklocalnum
911  * Array of the local indices of the cblk.
912  *
913  * @param[inout] solvcblk
914  * Pointer to the current cblk to register.
915  *
916  * @param[inout] solvblok
917  * Pointer to the first block of the current cblk.
918  *
919  * @param[in] lcblknm
920  * The local index of the cblk.
921  *
922  * @param[in] brownum
923  * The current index in the browtab.
924  *
925  * @param[in] gcblknm
926  * The global index of the current cblk.
927  *
928  * @param[in] ownerid
929  * The index of the local MPI rank.
930  *
931  * @return The pointer to the next solver block to register.
932  *
933  *******************************************************************************/
934 SolverBlok*
936  const Cand *candcblk,
937  const pastix_int_t *cblklocalnum,
938  SolverCblk *solvcblk,
939  SolverBlok *solvblok,
940  pastix_int_t lcblknm,
941  pastix_int_t brownum,
942  pastix_int_t gcblknm,
943  pastix_int_t ownerid )
944 {
945  symbol_cblk_t *symbcblk = symbmtx->cblktab + gcblknm;
946  symbol_blok_t *symbblok = symbmtx->bloktab + symbcblk->bloknum;
947  SolverBlok *fblokptr = solvblok;
948  pastix_int_t stride = 0;
949  pastix_int_t layout2D = candcblk->cblktype & CBLK_LAYOUT_2D;
950  pastix_int_t fcolnum, lcolnum, nbcols, j;
951 
952  assert( solvblok != NULL );
953  assert( brownum >= 0 );
954  assert( symbblok->lcblknm == gcblknm );
955  assert( (cblklocalnum == NULL) || (lcblknm == cblklocalnum[gcblknm]) );
956 
957  /*
958  * Compute the number of columns of the fan-in
959  */
960  nbcols = symbol_cblk_get_colnum( symbmtx, symbcblk, &fcolnum, &lcolnum );
961 
962  /*
963  * Register all the local blocks
964  */
965  for ( j=symbcblk[0].bloknum; j<symbcblk[1].bloknum; j++, symbblok++ )
966  {
967  pastix_int_t frownum, lrownum, nbrows;
968 
969  nbrows = symbol_blok_get_rownum( symbmtx, symbblok, &frownum, &lrownum );
970  assert( nbrows >= 1 );
971 
972  /* Init the blok */
973  solvMatGen_init_blok( solvblok, lcblknm,
974  cblklocalnum == NULL ? symbblok->fcblknm : cblklocalnum[symbblok->fcblknm],
975  frownum, lrownum, stride, nbcols,
976  layout2D );
977  solvblok->gbloknm = j;
978  stride += nbrows;
979  solvblok++;
980  }
981 
982  solvMatGen_init_cblk( solvcblk, fblokptr, candcblk, symbcblk,
983  fcolnum, lcolnum, brownum, stride,
984  gcblknm, ownerid );
985 
986  solvcblk->lcolidx = fcolnum;
987 
988  return solvblok;
989 }
990 
991 /**
992  *******************************************************************************
993  *
994  * @brief Register a remote cblk from a solver_recv_cblk_t structure (Fanin|Recv)
995  *
996  *******************************************************************************
997  *
998  * @param[in] symbmtx
999  * The pointer to the symbol matrix.
1000  *
1001  * @param[in] recvcblk
1002  * The associated solver_recv_cblk_t structure used to initialize the
1003  * remote cblk.
1004  *
1005  * @param[in] candcblk
1006  * The cand structure associated to the current cblk to get the type of
1007  * the cblk.
1008  *
1009  * @param[in] cblklocalnum
1010  * Array of the local indices of the cblk.
1011  *
1012  * @param[inout] solvcblk
1013  * Pointer to the current cblk to register.
1014  *
1015  * @param[inout] solvblok
1016  * Pointer to the first block of the current cblk.
1017  *
1018  * @param[in] lcblknm
1019  * The local index of the cblk.
1020  *
1021  * @param[in] brownum
1022  * The current index in the browtab.
1023  *
1024  * @param[in] gcblknm
1025  * The global index of the current cblk.
1026  *
1027  * @return The pointer to the next solver block to register.
1028  *
1029  *******************************************************************************/
1030 SolverBlok*
1032  const solver_cblk_recv_t *recvcblk,
1033  const Cand *candcblk,
1034  const pastix_int_t *cblklocalnum,
1035  SolverCblk *solvcblk,
1036  SolverBlok *solvblok,
1037  pastix_int_t lcblknm,
1038  pastix_int_t brownum,
1039  pastix_int_t gcblknm )
1040 {
1041  const solver_blok_recv_t *recvblok = recvcblk->bloktab;
1042  symbol_cblk_t *symbcblk = symbmtx->cblktab + gcblknm;
1043  symbol_blok_t *symbblok = symbmtx->bloktab + symbcblk->bloknum;
1044  SolverBlok *fblokptr = solvblok;
1045  pastix_int_t stride = 0;
1046  pastix_int_t layout2D = candcblk->cblktype & CBLK_LAYOUT_2D;
1047  pastix_int_t fcolnum, lcolnum, nbcols, j;
1048 
1049  assert( solvblok != NULL );
1050  assert( brownum >= 0 );
1051  assert( symbblok->lcblknm == gcblknm );
1052 
1053  /*
1054  * Compute the number of columns of the fan-in
1055  */
1056  if ( symbmtx->dof < 0 ) {
1057  fcolnum = symbmtx->dofs[recvcblk->fcolnum];
1058  lcolnum = symbmtx->dofs[recvcblk->lcolnum + 1] - 1;
1059  }
1060  else {
1061  fcolnum = symbmtx->dof * recvcblk->fcolnum;
1062  lcolnum = symbmtx->dof * ( recvcblk->lcolnum + 1 ) - 1;
1063  }
1064  nbcols = lcolnum - fcolnum + 1;
1065 
1066  /*
1067  * Register all the local blocks
1068  */
1069  for ( j=symbcblk[0].bloknum; j<symbcblk[1].bloknum; j++, recvblok++, symbblok++ )
1070  {
1071  pastix_int_t frownum, lrownum, nbrows;
1072 
1073  if ( symbmtx->dof < 0 ) {
1074  frownum = symbmtx->dofs[recvblok->frownum];
1075  lrownum = symbmtx->dofs[recvblok->lrownum + 1] - 1;
1076  }
1077  else {
1078  frownum = symbmtx->dof * recvblok->frownum;
1079  lrownum = symbmtx->dof * ( recvblok->lrownum + 1 ) - 1;
1080  }
1081  nbrows = lrownum - frownum + 1;
1082 
1083  if ( nbrows < 1 ) {
1084  continue;
1085  }
1086 
1087  /* Init the blok */
1088  solvMatGen_init_blok( solvblok,
1089  lcblknm, cblklocalnum[symbblok->fcblknm],
1090  frownum, lrownum, stride, nbcols,
1091  layout2D );
1092  solvblok->gbloknm = j;
1093  stride += nbrows;
1094  solvblok++;
1095  }
1096 
1097  solvMatGen_init_cblk( solvcblk, fblokptr, candcblk, symbcblk,
1098  fcolnum, lcolnum, brownum, stride,
1099  gcblknm, recvcblk->ownerid );
1100 
1101  solvcblk->lcolidx = -1;
1102 
1103 #if defined(PASTIX_BLEND_FANIN_FR)
1104  if( solvcblk->cblktype & CBLK_COMPRESSED ) {
1105  solvcblk->cblktype &= (~CBLK_COMPRESSED);
1106  }
1107 #endif
1108  /* No Schur complement in distributed for the moment */
1109  if( solvcblk->cblktype & CBLK_IN_SCHUR ) {
1110  solvcblk->cblktype &= (~CBLK_IN_SCHUR);
1111  }
1112 
1113  return solvblok;
1114 }
1115 
1116 /**
1117  *@}
1118  */
BEGIN_C_DECLS typedef int pastix_int_t
Definition: datatypes.h:51
int8_t cblktype
Definition: cand.h:36
Processor candidate group to own a column blok.
Definition: cand.h:28
pastix_int_t extendint_Size(const ExtendVectorINT *)
Return the number of element stored in the vector.
Definition: extendVector.c:112
pastix_int_t extendint_Read(const ExtendVectorINT *, pastix_int_t)
Return the element of index eltnum.
Definition: extendVector.c:136
SimuCluster * clustab
Definition: simu.h:123
SimuProc * proctab
Definition: simu.h:122
int8_t owned
Definition: simu.h:80
SimuTask * tasktab
Definition: simu.h:121
int ownerclust
Definition: simu.h:93
ExtendVectorINT * tasktab
Definition: simu.h:60
pastix_int_t ctrbcnt
Definition: simu.h:108
pastix_int_t cblknum
Definition: simu.h:101
SimuBlok * bloktab
Definition: simu.h:126
pastix_int_t fprocnum
Definition: simu.h:47
pastix_int_t bloknum
Definition: simu.h:102
SimuCblk * cblktab
Definition: simu.h:125
pastix_int_t prionum
Definition: simu.h:100
pastix_int_t tasknbr
Definition: simu.h:119
Thread structure for the simulation.
Definition: simu.h:56
Task structure for the simulation.
Definition: simu.h:99
Control structure for the simulation.
Definition: simu.h:116
SolverBlok * solvMatGen_register_local_cblk(const symbol_matrix_t *symbmtx, const Cand *candcblk, const pastix_int_t *cblklocalnum, SolverCblk *solvcblk, SolverBlok *solvblok, pastix_int_t lcblknm, pastix_int_t brownum, pastix_int_t gcblknm, pastix_int_t ownerid)
Register a local cblk from a symbol_cblk_t structure !(Fanin|Recv)
void solvMatGen_fill_localnums(const symbol_matrix_t *symbmtx, const SimuCtrl *simuctrl, SolverMatrix *solvmtx, pastix_int_t *cblklocalnum, pastix_int_t *bloklocalnum, pastix_int_t *tasklocalnum, solver_cblk_recv_t **ftgttab)
Fill the local numbering arrays to compress the symbol information into solver.
int solver_recv_get_bloknbr(const solver_cblk_recv_t *ftgtptr, const symbol_cblk_t *symbcblk, const symbol_blok_t *symbblok)
Compute the number of valid blocks in fanin/recv cblk.
Definition: solver_recv.c:306
void solver_recv_update_fanin(solver_cblk_recv_t **faninptr, const symbol_matrix_t *symbmtx, const symbol_cblk_t *cblk, const symbol_blok_t *blok, const symbol_cblk_t *fcblk, int ownerid)
Register a new contribution to a fanin cblk.
Definition: solver_recv.c:203
SolverBlok * solvMatGen_register_remote_cblk(const symbol_matrix_t *symbmtx, const solver_cblk_recv_t *recvcblk, const Cand *candcblk, const pastix_int_t *cblklocalnum, SolverCblk *solvcblk, SolverBlok *solvblok, pastix_int_t lcblknm, pastix_int_t brownum, pastix_int_t gcblknm)
Register a remote cblk from a solver_recv_cblk_t structure (Fanin|Recv)
void solver_recv_update_recv(solver_cblk_recv_t **recvptr, const symbol_matrix_t *symbmtx, const symbol_cblk_t *cblk, const symbol_blok_t *blok, const symbol_cblk_t *fcblk, int ownerid)
Register a new contribution to a recv cblk.
Definition: solver_recv.c:250
static void solvMatGen_init_cblk(SolverCblk *solvcblk, SolverBlok *fblokptr, const Cand *candcblk, const symbol_cblk_t *symbcblk, pastix_int_t fcolnum, pastix_int_t lcolnum, pastix_int_t brownum, pastix_int_t stride, pastix_int_t cblknum, int ownerid)
Initialize a solver cblk.
void solvMatGen_fill_tasktab(SolverMatrix *solvmtx, isched_t *isched, const SimuCtrl *simuctrl, const pastix_int_t *tasklocalnum, const pastix_int_t *cblklocalnum, const pastix_int_t *bloklocalnum, pastix_int_t clustnum, int is_dbg)
Fill the global tasktab array, as well as the thread ttsktab arrays.
void solvMatGen_stats_last(SolverMatrix *solvmtx)
Mark blocks if they belong to the last supernode, or if they are facing it for statistical purpose on...
void solvMatGen_fill_ttsktab(isched_thread_t *ctx, void *args)
Fill the ttsktab for it's own thread.
void solvMatGen_fill_ttsktab_dbg(isched_thread_t *ctx, void *args)
Fill in ttsktab for it's own thread. Only for debugging factorization.
pastix_int_t solvMatGen_reorder_browtab(const symbol_matrix_t *symbmtx, const symbol_cblk_t *symbcblk, SolverMatrix *solvmtx, SolverCblk *solvcblk, pastix_int_t *browtmp, const pastix_int_t *cblklocalnum, const pastix_int_t *bloklocalnum, pastix_int_t brownum)
Reorder the browtab from the symbol structure in a distributed way. First stock the 1D blocks and the...
void solvMatGen_max_buffers(SolverMatrix *solvmtx)
Compute the maximum area of the temporary buffers used during computation.
static void solvMatGen_init_blok(SolverBlok *solvblok, pastix_int_t lcblknm, pastix_int_t fcblknm, pastix_int_t frownum, pastix_int_t lrownum, pastix_int_t stride, pastix_int_t nbcols, pastix_int_t layout2D)
Initialize a solver block.
pastix_int_t fcblknm
Definition: symbol.h:63
pastix_int_t frownum
Definition: symbol.h:60
pastix_int_t lrownum
Definition: symbol.h:61
pastix_int_t * browtab
Definition: symbol.h:85
symbol_cblk_t * cblktab
Definition: symbol.h:83
pastix_int_t bloknum
Definition: symbol.h:48
pastix_int_t brownum
Definition: symbol.h:49
pastix_int_t dof
Definition: symbol.h:87
symbol_blok_t * bloktab
Definition: symbol.h:84
pastix_int_t lcblknm
Definition: symbol.h:62
pastix_int_t * dofs
Definition: symbol.h:89
pastix_int_t cblknbr
Definition: symbol.h:79
static pastix_int_t symbol_blok_get_rownum(const symbol_matrix_t *symbmtx, symbol_blok_t *symbblok, pastix_int_t *frownum, pastix_int_t *lrownum)
Get the expanded row index of a symbol_blok.
Definition: symbol.h:155
static pastix_int_t symbol_cblk_get_colnum(const symbol_matrix_t *symbmtx, symbol_cblk_t *symbcblk, pastix_int_t *fcolnum, pastix_int_t *lcolnum)
Get the expanded column indexes of a symbol_cblk.
Definition: symbol.h:116
Symbol block structure.
Definition: symbol.h:59
Symbol column block structure.
Definition: symbol.h:45
Symbol matrix structure.
Definition: symbol.h:77
pastix_int_t lrownum
Definition: solver.h:99
pastix_int_t browind
Definition: solver.h:145
static pastix_int_t blok_rownbr(const SolverBlok *blok)
Compute the number of rows of a block.
Definition: solver.h:390
pastix_int_t taskid
Definition: solver.h:119
pastix_int_t cblknum
Definition: solver.h:121
pastix_int_t lrownum
Definition: solver.h:143
pastix_int_t brownum
Definition: solver.h:166
pastix_int_t priority
Definition: solver.h:177
static pastix_int_t cblk_colnbr(const SolverCblk *cblk)
Compute the number of columns in a column block.
Definition: solver.h:324
pastix_int_t brown2d
Definition: solver.h:167
pastix_int_t fcblknm
Definition: solver.h:140
pastix_int_t frownum
Definition: solver.h:98
pastix_int_t cblknbr
Definition: solver.h:208
pastix_int_t fcolnum
Definition: solver.h:110
SolverBlok *restrict bloktab
Definition: solver.h:223
pastix_int_t frownum
Definition: solver.h:142
pastix_int_t brownbr
Definition: solver.h:221
pastix_int_t gbloknm
Definition: solver.h:141
pastix_int_t lcolnum
Definition: solver.h:111
pastix_int_t faninnbr
Definition: solver.h:209
pastix_int_t prionum
Definition: solver.h:120
SolverBlok * fblokptr
Definition: solver.h:163
pastix_int_t bloknbr
Definition: solver.h:220
pastix_int_t recvnbr
Definition: solver.h:212
pastix_int_t *restrict browtab
Definition: solver.h:224
pastix_int_t volatile ctrbcnt
Definition: solver.h:123
pastix_int_t lcblknm
Definition: solver.h:139
solver_blok_recv_t bloktab[1]
Definition: solver.h:112
int threadid
Definition: solver.h:176
pastix_int_t lcolidx
Definition: solver.h:165
int8_t inlast
Definition: solver.h:146
SolverCblk *restrict cblktab
Definition: solver.h:222
pastix_int_t stride
Definition: solver.h:164
int8_t cblktype
Definition: solver.h:159
int ownerid
Definition: solver.h:175
pastix_int_t bloknum
Definition: solver.h:122
Solver recv block structure.
Definition: solver.h:97
Solver block structure.
Definition: solver.h:137
Solver recv column block structure.
Definition: solver.h:107
Solver column block structure.
Definition: solver.h:156
Solver column block structure.
Definition: solver.h:200
The task structure for the numerical factorization.
Definition: solver.h:118