PaStiX Handbook  6.2.1
pastix_task_sopalin.c
Go to the documentation of this file.
1 /**
2  *
3  * @file pastix_task_sopalin.c
4  *
5  * PaStiX factorization routines
6  *
7  * @copyright 2004-2021 Bordeaux INP, CNRS (LaBRI UMR 5800), Inria,
8  * Univ. Bordeaux. All rights reserved.
9  *
10  * @version 6.2.1
11  * @author Pascal Henon
12  * @author Xavier Lacoste
13  * @author Pierre Ramet
14  * @author Mathieu Faverge
15  * @author Esragul Korkmaz
16  * @author Gregoire Pichon
17  * @author Theophile Terraz
18  * @date 2021-07-02
19  *
20  **/
21 #include "common.h"
22 #include <lapacke.h>
23 #include "isched.h"
24 #include <spm.h>
25 #include "bcsc/bcsc.h"
26 #include "blend/solver.h"
27 #include "coeftab.h"
28 #include "sopalin/sopalin_data.h"
29 #include "kernels/pastix_lowrank.h"
34 #include "kernels/kernels_trace.h"
35 
36 #if defined(PASTIX_WITH_PARSEC)
38 #endif
39 
40 #if defined(PASTIX_WITH_STARPU)
42 #endif
43 
44 static void (*sopalinFacto[5][4])(pastix_data_t *, sopalin_data_t*) =
45 {
46  { sopalin_spotrf, sopalin_dpotrf, sopalin_cpotrf, sopalin_zpotrf },
47  { sopalin_ssytrf, sopalin_dsytrf, sopalin_csytrf, sopalin_zsytrf },
48  { sopalin_sgetrf, sopalin_dgetrf, sopalin_cgetrf, sopalin_zgetrf },
49  { sopalin_spotrf, sopalin_dpotrf, sopalin_cpxtrf, sopalin_zpxtrf },
50  { sopalin_ssytrf, sopalin_dsytrf, sopalin_chetrf, sopalin_zhetrf }
51 };
52 
53 /**
54  *******************************************************************************
55  *
56  * @ingroup pastix_numfact
57  *
58  * @brief Fill the internal block CSC structure.
59  *
60  * This internal block CSC can be used by the refinement step with or without
61  * the preconditioner obtained by the numerical factorization.
62  *
63  * This routine is affected by the following parameters:
64  * IPARM_VERBOSE.
65  *
66  *******************************************************************************
67  *
68  * @param[inout] pastix_data
69  * The pastix_data structure that describes the solver instance.
70  * On exit, the internal block CSC is filled with entries from
71  * the spm matrix.
72  *
73  * @param[inout] spm
74  * The sparse matrix descriptor that describes problem instance.
75  * On exit, if IPARM_FREE_CSCUSER is set, the spm is freed.
76  *
77  *******************************************************************************
78  *
79  * @retval PASTIX_SUCCESS on successful exit,
80  * @retval PASTIX_ERR_BADPARAMETER if one parameter is incorrect,
81  * @retval PASTIX_ERR_OUTOFMEMORY if one allocation failed.
82  *
83  *******************************************************************************/
84 int
85 pastix_subtask_spm2bcsc( pastix_data_t *pastix_data,
86  spmatrix_t *spm )
87 {
88  double time;
89 
90  /*
91  * Check parameters
92  */
93  if (pastix_data == NULL) {
94  errorPrint("pastix_subtask_spm2bcsc: wrong pastix_data parameter");
96  }
97  if (spm == NULL) {
98  errorPrint("pastix_subtask_spm2bcsc: wrong spm parameter");
100  }
101  if ( !(pastix_data->steps & STEP_ANALYSE) ) {
102  errorPrint("pastix_subtask_spm2bcsc: All steps from pastix_task_init() to pastix_task_blend() have to be called before calling this function");
104  }
105 
106  /*
107  * Compute the norm of A, to scale the epsilon parameter for pivoting
108  */
109  {
110  pastix_data->dparm[ DPARM_A_NORM ] = spmNorm( SpmFrobeniusNorm, spm );
111  if ( pastix_data->iparm[IPARM_VERBOSE] > PastixVerboseNo ) {
112  pastix_print( pastix_data->inter_node_procnum, 0,
113  " ||A||_2 = %e\n",
114  pastix_data->dparm[ DPARM_A_NORM ] );
115  }
116  }
117 
118  /*
119  * Switch the SolverMatrix pointer according to the runtime.
120  */
121  {
122  pastix_int_t isched = pastix_data->iparm[IPARM_SCHEDULER];
123 
124  if ( isSchedRuntime( isched ) ) {
125  pastix_data->solvmatr = pastix_data->solvglob;
126  }
127  else {
128  assert( isSchedPthread( isched ) );
129  pastix_data->solvmatr = pastix_data->solvloc;
130  }
131  pastix_data->sched = isched;
132  }
133 
134  /*
135  * Fill in the internal blocked CSC. We consider that if this step is called
136  * the spm values have changed so we need to update the blocked csc.
137  */
138  if ( pastix_data->bcsc != NULL )
139  {
140  bcscExit( pastix_data->bcsc );
141  memFree_null( pastix_data->bcsc );
142  }
143 
144  MALLOC_INTERN( pastix_data->bcsc, 1, pastix_bcsc_t );
145 
146  time = bcscInit( spm,
147  pastix_data->ordemesh,
148  pastix_data->solvmatr,
149  (pastix_data->iparm[IPARM_FACTORIZATION] == PastixFactLU), /*&& (! pastix_data->iparm[IPARM_ONLY_REFINE]) )*/
150  pastix_data->bcsc );
151 
152  if ( pastix_data->iparm[IPARM_VERBOSE] > PastixVerboseNot ) {
153  pastix_print( pastix_data->inter_node_procnum, 0, OUT_BCSC_TIME, time );
154  }
155 
156  if ( pastix_data->iparm[IPARM_FREE_CSCUSER] ) {
157  spmExit( spm );
158  }
159 
160  /*
161  * Invalidate following step, and add current step to the ones performed
162  */
163  pastix_data->steps &= ~STEP_BCSC2CTAB;
164  pastix_data->steps |= STEP_CSC2BCSC;
165 
166  return PASTIX_SUCCESS;
167 }
168 
169 /**
170  *******************************************************************************
171  *
172  * @ingroup pastix_numfact
173  *
174  * @brief Fill the internal solver matrix structure.
175  *
176  * This step is linked with the pastix_subtask_sopalin() since this structure is
177  * only used during the numerical factorization.
178  *
179  * This routine is affected by the following parameters:
180  * IPARM_VERBOSE, [IPARM_FACTORIZATION.
181  *
182  *******************************************************************************
183  *
184  * @param[inout] pastix_data
185  * The pastix_data structure that describes the solver instance.
186  * On exit, the internal solver structure is filled with entries from
187  * the internal block CSC matrix.
188  *
189  *******************************************************************************
190  *
191  * @retval PASTIX_SUCCESS on successful exit,
192  * @retval PASTIX_ERR_BADPARAMETER if one parameter is incorrect,
193  * @retval PASTIX_ERR_OUTOFMEMORY if one allocation failed.
194  *
195  *******************************************************************************/
196 int
197 pastix_subtask_bcsc2ctab( pastix_data_t *pastix_data )
198 {
199  pastix_bcsc_t *bcsc;
200  pastix_lr_t *lr;
201  Clock timer;
202  int mtxtype;
203 
204  /*
205  * Check parameters
206  */
207  if (pastix_data == NULL) {
208  errorPrint("pastix_subtask_bcsc2ctab: wrong pastix_data parameter");
210  }
211  if ( !(pastix_data->steps & STEP_CSC2BCSC) ) {
212  errorPrint("pastix_subtask_bcsc2ctab: All steps from pastix_task_init() to pastix_stask_blend() have to be called before calling this function");
214  }
215  if (pastix_data->bcsc == NULL) {
216  errorPrint("pastix_subtask_bcsc2ctab: wrong pastix_data->bcsc parameter");
218  }
219 
220  /* Ensure that the scheduler is correct */
221  pastix_check_and_correct_scheduler( pastix_data );
222 
223  clockSyncStart( timer, pastix_data->inter_node_comm );
224 
225  /* Initialize low-rank parameters */
226  lr = &(pastix_data->solvmatr->lowrank);
227  lr->compress_when = pastix_data->iparm[IPARM_COMPRESS_WHEN];
228  lr->compress_method = pastix_data->iparm[IPARM_COMPRESS_METHOD];
229  lr->compress_min_width = pastix_data->iparm[IPARM_COMPRESS_MIN_WIDTH];
230  lr->compress_min_height = pastix_data->iparm[IPARM_COMPRESS_MIN_HEIGHT];
231  lr->compress_preselect = pastix_data->iparm[IPARM_COMPRESS_PRESELECT];
232  lr->use_reltol = pastix_data->iparm[IPARM_COMPRESS_RELTOL];
233  lr->tolerance = pastix_data->dparm[DPARM_COMPRESS_TOLERANCE];
234  lr->ilu_lvl = pastix_data->iparm[IPARM_COMPRESS_ILUK];
235 
236  pastix_lr_minratio = pastix_data->dparm[DPARM_COMPRESS_MIN_RATIO];
237  pastix_lr_ortho = pastix_data->iparm[IPARM_COMPRESS_ORTHO];
238 
239  bcsc = pastix_data->bcsc;
240  lr->core_ge2lr = ge2lrMethods[ pastix_data->iparm[IPARM_COMPRESS_METHOD] ][bcsc->flttype-2];
241  lr->core_rradd = rraddMethods[ pastix_data->iparm[IPARM_COMPRESS_METHOD] ][bcsc->flttype-2];
242 
243  /*
244  * Update the ilu level based on when
245  */
246  if ( lr->ilu_lvl == -2 ) {
247  if ( pastix_data->iparm[IPARM_COMPRESS_WHEN] == PastixCompressWhenBegin ) {
248  lr->ilu_lvl = -1;
249  }
250  else if ( pastix_data->iparm[IPARM_COMPRESS_WHEN] == PastixCompressWhenEnd ) {
251  lr->ilu_lvl = INT_MAX;
252  }
253  else if ( pastix_data->iparm[IPARM_COMPRESS_WHEN] == PastixCompressWhenDuring ) {
254  lr->ilu_lvl = 0;
255  }
256  }
257 
258  /*
259  * Set the maximum rank function
260  */
261  if ( pastix_data->iparm[IPARM_COMPRESS_WHEN] == PastixCompressWhenBegin ) {
263  }
264  else {
266  }
267 
268  pastix_data->solvmatr->factotype = pastix_data->iparm[IPARM_FACTORIZATION];
269 
270 #if defined(PASTIX_WITH_MPI)
271  if ( ( pastix_data->iparm[IPARM_COMPRESS_WHEN] != PastixCompressNever ) &&
272  ( ( pastix_data->iparm[IPARM_SCHEDULER] == PastixSchedParsec ) ||
273  ( pastix_data->iparm[IPARM_SCHEDULER] == PastixSchedStarPU ) ) )
274  {
275  errorPrint( "pastix_task_sopalin: Low-Rank with MPI communication is not available yet with runtime" );
277  }
278 #endif
279  /*
280  * Fill in the internal coeftab structure. We consider that if this step is
281  * called the bcsc values have changed, or a factorization have already been
282  * performed, so we need to update the coeftab arrays.
283  */
284  if (pastix_data->bcsc != NULL)
285  {
286  coeftabExit( pastix_data->solvmatr );
287  }
288 
289  coeftabInit( pastix_data,
290  pastix_data->iparm[IPARM_FACTORIZATION] == PastixFactLU ? PastixLUCoef : PastixLCoef );
291 
292  switch( pastix_data->iparm[IPARM_FACTORIZATION] ) {
293  case PastixFactLLH:
294  case PastixFactLDLH:
295  mtxtype = PastixHermitian;
296  break;
297 
298  case PastixFactLLT:
299  case PastixFactLDLT:
300  mtxtype = PastixHermitian;
301  break;
302 
303  case PastixFactLU:
304  default:
305  mtxtype = PastixGeneral;
306  }
307 
308 #if defined(PASTIX_WITH_PARSEC)
309  if ( pastix_data->iparm[IPARM_SCHEDULER] == PastixSchedParsec )
310  {
311  /* Start PaRSEC if not already started */
312  if (pastix_data->parsec == NULL) {
313  int argc = 0;
314  pastix_parsec_init( pastix_data, &argc, NULL, NULL );
315  }
316  /* Create the matrix descriptor */
317  parsec_sparse_matrix_init( pastix_data->solvmatr,
318  pastix_size_of( bcsc->flttype ), mtxtype,
319  pastix_data->inter_node_procnbr,
320  pastix_data->inter_node_procnum );
321  }
322 #endif
323 
324 #if defined(PASTIX_WITH_STARPU)
325  if ( pastix_data->iparm[IPARM_SCHEDULER] == PastixSchedStarPU )
326  {
327  /* Create the matrix descriptor */
328  starpu_sparse_matrix_init( pastix_data->solvmatr,
329  pastix_size_of( bcsc->flttype ), mtxtype,
330  pastix_data->inter_node_procnbr,
331  pastix_data->inter_node_procnum );
332  }
333 #endif
334 
335  clockSyncStop( timer, pastix_data->inter_node_comm );
336  if ( pastix_data->iparm[IPARM_VERBOSE] > PastixVerboseNot ) {
337  pastix_print( pastix_data->inter_node_procnum, 0, OUT_COEFTAB_TIME,
338  clockVal(timer) );
339  }
340 
341  /* Invalidate following step, and add current step to the ones performed */
342  pastix_data->steps &= ~STEP_NUMFACT;
343  pastix_data->steps |= STEP_BCSC2CTAB;
344 
345  (void)mtxtype;
346  return PASTIX_SUCCESS;
347 }
348 
349 /**
350  *******************************************************************************
351  *
352  * @ingroup pastix_numfact
353  *
354  * @brief Factorize the given problem using Cholesky or LU decomposition.
355  *
356  * The user can call the pastix_task_solve() to obtain the solution.
357  *
358  * This routine is affected by the following parameters:
359  * IPARM_VERBOSE, IPARM_FACTORIZATION, IPARM_COMPRESS_WHEN.
360  *
361  *******************************************************************************
362  *
363  * @param[inout] pastix_data
364  * The pastix_data structure that describes the solver instance.
365  * On exit, the solver matrix structure stores the factorization of the
366  * given problem.
367  *
368  *******************************************************************************
369  *
370  * @retval PASTIX_SUCCESS on successful exit,
371  * @retval PASTIX_ERR_BADPARAMETER if one parameter is incorrect,
372  * @retval PASTIX_ERR_OUTOFMEMORY if one allocation failed.
373  *
374  *******************************************************************************/
375 int
376 pastix_subtask_sopalin( pastix_data_t *pastix_data )
377 {
378  sopalin_data_t sopalin_data;
379  SolverBackup_t *sbackup;
380  pastix_bcsc_t *bcsc;
381  PASTIX_Comm pastix_comm;
382  pastix_int_t *iparm;
383  double *dparm;
384 
385  /*
386  * Check parameters
387  */
388  if (pastix_data == NULL) {
389  errorPrint("pastix_subtask_sopalin: wrong pastix_data parameter");
391  }
392  if ( !(pastix_data->steps & STEP_ANALYSE) ) {
393  errorPrint("pastix_subtask_sopalin: All steps from pastix_task_init() to pastix_task_analyze() have to be called before calling this function");
395  }
396  if ( !(pastix_data->steps & STEP_CSC2BCSC) ) {
397  errorPrint("pastix_subtask_sopalin: All steps from pastix_task_init() to pastix_task_analyze() have to be called before calling this function");
399  }
400  if ( !(pastix_data->steps & STEP_BCSC2CTAB) ) {
401  errorPrint("pastix_subtask_sopalin: All steps from pastix_task_init() to pastix_task_analyze() have to be called before calling this function");
403  }
404 
405  bcsc = pastix_data->bcsc;
406  if (bcsc == NULL) {
407  errorPrint("pastix_subtask_sopalin: wrong pastix_data_bcsc parameter");
409  }
410 
411  /* Ensure that the scheduler is correct */
412  pastix_check_and_correct_scheduler( pastix_data );
413 
414  pastix_comm = pastix_data->inter_node_comm;
415  iparm = pastix_data->iparm;
416  dparm = pastix_data->dparm;
417 
418  /* Prepare the sopalin_data structure */
419  {
420  double threshold;
421 
422  sopalin_data.solvmtx = pastix_data->solvmatr;
423 
424  /* TODO: might change the behavior: if the user wants a ratio of the norm, it could compute it himself */
425  if ( dparm[ DPARM_EPSILON_MAGN_CTRL ] < 0. ) {
426  threshold = - dparm[ DPARM_EPSILON_MAGN_CTRL ];
427  }
428  else if ( dparm[ DPARM_EPSILON_MAGN_CTRL ] == 0. ) {
429  /*
430  * Use the rule presented in "Making Sparse Gaussian Elimination
431  * Scalable by Static Pivoting", S. X. Li, J. Demmel
432  *
433  * sqrt(eps) * ||A||_1 is a small half precision pertrbation
434  * intriduced in the pb when the pivot is too small
435  */
436  double eps;
437  if ( (bcsc->flttype == PastixFloat) || (bcsc->flttype == PastixComplex32) ) {
438  eps = LAPACKE_slamch_work( 'e' );
439  }
440  else {
441  eps = LAPACKE_dlamch_work( 'e' );
442  }
443  threshold = sqrt(eps) * dparm[DPARM_A_NORM];
444  }
445  else {
446  threshold = dparm[ DPARM_EPSILON_MAGN_CTRL ] * dparm[DPARM_A_NORM];
447  }
448 
449  sopalin_data.solvmtx->diagthreshold = threshold;
450  sopalin_data.solvmtx->nbpivots = 0;
451 
452  sopalin_data.cpu_coefs = &(pastix_data->cpu_models->coefficients[bcsc->flttype-2]);
453  sopalin_data.gpu_coefs = &(pastix_data->gpu_models->coefficients[bcsc->flttype-2]);
454  }
455 
456  sbackup = solverBackupInit( pastix_data->solvmatr );
457  pastix_data->solvmatr->restore = 2;
458  {
459  void (*factofct)( pastix_data_t *, sopalin_data_t *);
460  double timer, timer_local, flops, flops_g = 0.;
461 
462  factofct = sopalinFacto[ pastix_data->iparm[IPARM_FACTORIZATION] ][bcsc->flttype-2];
463  assert(factofct);
464 
465  kernelsTraceStart( pastix_data );
466  clockSyncStart( timer, pastix_comm );
467  clockStart(timer_local);
468 
469  factofct( pastix_data, &sopalin_data );
470 
471  clockStop(timer_local);
472  clockSyncStop( timer, pastix_comm );
473  kernelsTraceStop( pastix_data );
474 
475  /* Output time and flops */
476  pastix_data->dparm[DPARM_FACT_TIME] = clockVal(timer);
477  flops = pastix_data->dparm[DPARM_FACT_THFLOPS] / pastix_data->dparm[DPARM_FACT_TIME];
478  pastix_data->dparm[DPARM_FACT_FLOPS] = ((flops / 1024.) / 1024.) / 1024.;
479 
480  /* Compute the number of flops performed (may be different from THFLOPS for low-rank) */
481  {
482  double flops_l = pastix_data->dparm[DPARM_FACT_RLFLOPS];
483  MPI_Reduce( &flops_l, &flops_g, 1, MPI_DOUBLE, MPI_SUM, 0, pastix_data->inter_node_comm );
484  }
485 
486  /* Reduce the number of pivots */
487  {
488  int nbpivot_l = sopalin_data.solvmtx->nbpivots;
489  int nbpivot_g;
490  MPI_Allreduce( &nbpivot_l, &nbpivot_g, 1, MPI_INT, MPI_SUM, pastix_data->inter_node_comm );
491  pastix_data->iparm[IPARM_STATIC_PIVOTING] = nbpivot_g;
492  }
493 
494  if ( iparm[IPARM_VERBOSE] > PastixVerboseNot ) {
495  pastix_print( pastix_data->inter_node_procnum, 0, OUT_SOPALIN_TIME,
496  clockVal(timer),
497  pastix_print_value( flops ),
498  pastix_print_unit( flops ),
499  pastix_print_value( flops_g ),
500  pastix_print_unit( flops_g ),
501  (long)pastix_data->iparm[IPARM_STATIC_PIVOTING] );
502  }
503 
504 #if defined(PASTIX_WITH_PARSEC) && defined(PASTIX_DEBUG_PARSEC)
505  {
506  int i;
507  fprintf(stderr, "-- Check status of PaRSEC nbtasks counters\n" );
508  for (i=0; i<32; i++) {
509  if ( parsec_nbtasks_on_gpu[i] != 0 ) {
510  fprintf(stderr, "Device %d => %d\n",
511  i, parsec_nbtasks_on_gpu[i] );
512  }
513  parsec_nbtasks_on_gpu[i] = 0;
514  }
515  }
516 #endif /* defined(PASTIX_WITH_PARSEC) && defined(PASTIX_DEBUG_PARSEC) */
517  }
518  solverBackupRestore( pastix_data->solvmatr, sbackup );
519  solverBackupExit( sbackup );
520 
521 #if defined(PASTIX_NUMFACT_DUMP_SOLVER)
522  {
523  FILE *stream = NULL;
524 
525  pastix_gendirectories( pastix_data );
526  stream = pastix_fopenw( pastix_data->dir_local, "solver.eps", "w" );
527  if ( stream ) {
528  solverDraw( pastix_data->solvmatr,
529  stream, iparm[IPARM_VERBOSE],
530  pastix_data->dir_local );
531  fclose(stream);
532  }
533  }
534 #endif
535 
536  if ( (pastix_data->iparm[IPARM_VERBOSE] > PastixVerboseNot) &&
537  (pastix_data->iparm[IPARM_COMPRESS_WHEN] != PastixCompressNever) )
538  {
539  /* Compute the memory gain */
540  coeftabMemory[bcsc->flttype-2]( pastix_data->solvmatr );
541  }
542 
543  /* Invalidate following steps, and add factorization step to the ones performed */
544  pastix_data->steps &= ~( STEP_BCSC2CTAB |
545  STEP_SOLVE |
546  STEP_REFINE );
547  pastix_data->steps |= STEP_NUMFACT;
548 
549  (void)pastix_comm;
550  return EXIT_SUCCESS;
551 }
552 
553 /**
554  *******************************************************************************
555  *
556  * @ingroup pastix_users
557  *
558  * @brief Perform all the numerical factorization steps:
559  * fill the internal block CSC and the solver matrix structures,
560  * then apply the factorization step.
561  *
562  * The user can call the pastix_task_solve() to obtain the solution.
563  *
564  * This routine is affected by the following parameters:
565  * IPARM_VERBOSE.
566  *
567  *******************************************************************************
568  *
569  * @param[inout] pastix_data
570  * The pastix_data structure that describes the solver instance.
571  * On exit, the internal block CSC is filled with entries from the
572  * spm matrix and the solver matrix structure stores the factorization
573  * of the given problem.
574  *
575  * @param[inout] spm
576  * The sparse matrix descriptor that describes problem instance.
577  * On exit, the spm structure is freed if IPARM_FREE_CSCUSER is true.
578  *
579  *******************************************************************************
580  *
581  * @retval PASTIX_SUCCESS on successful exit,
582  * @retval PASTIX_ERR_BADPARAMETER if one parameter is incorrect,
583  * @retval PASTIX_ERR_OUTOFMEMORY if one allocation failed.
584  *
585  *******************************************************************************/
586 int
587 pastix_task_numfact( pastix_data_t *pastix_data,
588  spmatrix_t *spm )
589 {
590  pastix_int_t *iparm;
591  pastix_int_t procnum;
592  pastix_int_t steps;
593  int rc;
594 
595  /*
596  * Check parameters
597  */
598  if (pastix_data == NULL) {
599  errorPrint("pastix_task_sopalin: wrong pastix_data parameter");
601  }
602  if (spm == NULL) {
603  errorPrint("pastix_task_sopalin: wrong spm parameter");
605  }
606  if ( !(pastix_data->steps & STEP_ANALYSE) ) {
607  errorPrint("pastix_task_sopalin: All steps from pastix_task_init() to pastix_task_blend() have to be called before calling this function");
609  }
610 
611  iparm = pastix_data->iparm;
612  procnum = pastix_data->inter_node_procnum;
613 
614  if (iparm[IPARM_VERBOSE] > PastixVerboseNot) {
615  pastix_print( procnum, 0, OUT_STEP_SOPALIN,
616  pastixFactotypeStr( iparm[IPARM_FACTORIZATION] ) );
617  }
618 
619  /* Invalidate upcoming steps */
620  steps = ~(STEP_CSC2BCSC |
621  STEP_BCSC2CTAB |
622  STEP_NUMFACT |
623  STEP_SOLVE |
624  STEP_REFINE );
625  pastix_data->steps &= steps;
626 
627  rc = pastix_subtask_spm2bcsc( pastix_data, spm );
628  if (rc != PASTIX_SUCCESS) {
629  return rc;
630  }
631 
632  rc = pastix_subtask_bcsc2ctab( pastix_data );
633  if (rc != PASTIX_SUCCESS) {
634  return rc;
635  }
636 
637  rc = pastix_subtask_sopalin( pastix_data );
638  if (rc != PASTIX_SUCCESS) {
639  return rc;
640  }
641 
642  return EXIT_SUCCESS;
643 }
pastix_lr_minratio
double pastix_lr_minratio
Define the minimal ratio for which we accept to compress a matrix into a low-rank form or not.
Definition: lowrank.c:24
solver.h
solverBackupExit
void solverBackupExit(SolverBackup_t *b)
Free the solver backup data structure.
Definition: solver_backup.c:207
PastixFactLLT
@ PastixFactLLT
Definition: api.h:303
core_get_rklimit
pastix_int_t(* core_get_rklimit)(pastix_int_t, pastix_int_t)
Compute the maximal rank accepted for a given matrix size. The pointer is set according to the low-ra...
Definition: kernels_trace.c:30
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
DPARM_FACT_FLOPS
@ DPARM_FACT_FLOPS
Definition: api.h:164
core_get_rklimit_end
static pastix_int_t core_get_rklimit_end(pastix_int_t M, pastix_int_t N)
Compute the maximal rank accepted for a given matrix size for Just-In-Time strategy.
Definition: pastix_lowrank.h:84
DPARM_COMPRESS_TOLERANCE
@ DPARM_COMPRESS_TOLERANCE
Definition: api.h:173
pastix_lowrank.h
pastix_bcsc_s
Internal column block distributed CSC matrix.
Definition: bcsc.h:37
bcscExit
void bcscExit(pastix_bcsc_t *bcsc)
Cleanup the block csc structure but do not free the bcsc pointer.
Definition: bcsc.c:437
pastix_lr_s::use_reltol
int use_reltol
Definition: pastix_lowrank.h:153
pastix_lr_s::ilu_lvl
int ilu_lvl
Definition: pastix_lowrank.h:154
pastix_lr_s::compress_min_width
pastix_int_t compress_min_width
Definition: pastix_lowrank.h:150
IPARM_COMPRESS_RELTOL
@ IPARM_COMPRESS_RELTOL
Definition: api.h:131
pastix_gendirectories
void pastix_gendirectories(pastix_data_t *pastix_data)
Generate a unique temporary directory to store output files.
Definition: api.c:69
pastix_parsec.h
IPARM_STATIC_PIVOTING
@ IPARM_STATIC_PIVOTING
Definition: api.h:100
pastix_lr_ortho
pastix_int_t pastix_lr_ortho
Define the orthogonalization method.
Definition: lowrank.c:25
IPARM_COMPRESS_ORTHO
@ IPARM_COMPRESS_ORTHO
Definition: api.h:130
PastixCompressNever
@ PastixCompressNever
Definition: api.h:362
pastix_lr_s::compress_method
pastix_compress_method_t compress_method
Definition: pastix_lowrank.h:149
ge2lrMethods
const fct_ge2lr_t ge2lrMethods[PastixCompressMethodNbr][4]
Array of pointers to the multiple arithmetic and algorithmic variants of ge2lr.
Definition: lowrank.c:43
PastixGeneral
@ PastixGeneral
Definition: api.h:433
DPARM_A_NORM
@ DPARM_A_NORM
Definition: api.h:172
IPARM_COMPRESS_MIN_HEIGHT
@ IPARM_COMPRESS_MIN_HEIGHT
Definition: api.h:127
pastix_lr_s::compress_preselect
int compress_preselect
Definition: pastix_lowrank.h:152
pastix_subtask_spm2bcsc
int pastix_subtask_spm2bcsc(pastix_data_t *pastix_data, spmatrix_t *spm)
Fill the internal block CSC structure.
Definition: pastix_task_sopalin.c:85
pastix_lr_s::compress_min_height
pastix_int_t compress_min_height
Definition: pastix_lowrank.h:151
IPARM_COMPRESS_METHOD
@ IPARM_COMPRESS_METHOD
Definition: api.h:129
pastix_bcsc_s::flttype
int flttype
Definition: bcsc.h:41
kernelsTraceStart
void kernelsTraceStart(const pastix_data_t *pastix_data)
Start the trace module.
Definition: kernels_trace.c:66
IPARM_SCHEDULER
@ IPARM_SCHEDULER
Definition: api.h:116
DPARM_FACT_THFLOPS
@ DPARM_FACT_THFLOPS
Definition: api.h:165
pastix_zlrcores.h
bcsc.h
PastixSchedParsec
@ PastixSchedParsec
Definition: api.h:313
pastix_lr_s::core_ge2lr
fct_ge2lr_t core_ge2lr
Definition: pastix_lowrank.h:157
PastixCompressWhenBegin
@ PastixCompressWhenBegin
Definition: api.h:363
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
IPARM_COMPRESS_MIN_WIDTH
@ IPARM_COMPRESS_MIN_WIDTH
Definition: api.h:126
PastixFactLLH
@ PastixFactLLH
Definition: api.h:300
PastixVerboseNot
@ PastixVerboseNot
Definition: api.h:207
PASTIX_SUCCESS
@ PASTIX_SUCCESS
Definition: api.h:344
PastixVerboseNo
@ PastixVerboseNo
Definition: api.h:208
pastix_fopenw
FILE * pastix_fopenw(const char *dirname, const char *filename, const char *mode)
Open a file in the unique directory of the pastix instance.
Definition: api.c:232
IPARM_COMPRESS_WHEN
@ IPARM_COMPRESS_WHEN
Definition: api.h:128
PastixLCoef
@ PastixLCoef
Definition: api.h:454
pastix_slrcores.h
IPARM_FREE_CSCUSER
@ IPARM_FREE_CSCUSER
Definition: api.h:101
parsec_sparse_matrix_init
void parsec_sparse_matrix_init(SolverMatrix *solvmtx, int typesize, int mtxtype, int nodes, int myrank)
Generate the PaRSEC descriptor of the sparse matrix.
Definition: parsec_sparse_matrix.c:551
PastixHermitian
@ PastixHermitian
Definition: api.h:435
DPARM_FACT_TIME
@ DPARM_FACT_TIME
Definition: api.h:163
DPARM_COMPRESS_MIN_RATIO
@ DPARM_COMPRESS_MIN_RATIO
Definition: api.h:174
starpu_sparse_matrix_init
void starpu_sparse_matrix_init(SolverMatrix *solvmtx, int typesize, int mtxtype, int nodes, int myrank)
Generate the StarPU descriptor of the sparse matrix.
Definition: starpu_sparse_matrix.c:163
IPARM_COMPRESS_ILUK
@ IPARM_COMPRESS_ILUK
Definition: api.h:133
coeftabExit
void coeftabExit(SolverMatrix *solvmtx)
Free the solver matrix structure.
Definition: coeftab.c:182
PastixCompressWhenDuring
@ PastixCompressWhenDuring
Definition: api.h:365
PastixFactLDLT
@ PastixFactLDLT
Definition: api.h:301
pastix_lr_s::tolerance
double tolerance
Definition: pastix_lowrank.h:155
coeftabMemory
coeftab_fct_memory_t coeftabMemory[4]
List of functions to compute the memory gain in low-rank per precision.
Definition: coeftab.c:39
PastixLUCoef
@ PastixLUCoef
Definition: api.h:456
coeftab.h
PastixCompressWhenEnd
@ PastixCompressWhenEnd
Definition: api.h:364
pastix_lr_s::core_rradd
fct_rradd_t core_rradd
Definition: pastix_lowrank.h:156
pastix_starpu.h
core_get_rklimit_begin
static pastix_int_t core_get_rklimit_begin(pastix_int_t M, pastix_int_t N)
Compute the maximal rank accepted for a given matrix size for Minimal-Memory strategy.
Definition: pastix_lowrank.h:95
IPARM_VERBOSE
@ IPARM_VERBOSE
Definition: api.h:36
rraddMethods
const fct_rradd_t rraddMethods[PastixCompressMethodNbr][4]
Array of pointers to the multiple arithmetic and algorithmic variants of rradd.
Definition: lowrank.c:52
PastixFactLU
@ PastixFactLU
Definition: api.h:302
pastix_subtask_bcsc2ctab
int pastix_subtask_bcsc2ctab(pastix_data_t *pastix_data)
Fill the internal solver matrix structure.
Definition: pastix_task_sopalin.c:197
pastix_lr_s::compress_when
pastix_compress_when_t compress_when
Definition: pastix_lowrank.h:148
solverBackupInit
SolverBackup_t * solverBackupInit(const SolverMatrix *solvmtx)
Initialize the backup structure.
Definition: solver_backup.c:59
coeftabInit
void coeftabInit(pastix_data_t *pastix_data, pastix_coefside_t side)
Initialize the solver matrix structure.
Definition: coeftab.c:135
pastix_dlrcores.h
DPARM_FACT_RLFLOPS
@ DPARM_FACT_RLFLOPS
Definition: api.h:166
bcscInit
double bcscInit(const spmatrix_t *spm, const pastix_order_t *ord, const SolverMatrix *solvmtx, pastix_int_t initAt, pastix_bcsc_t *bcsc)
Initialize the block csc matrix.
Definition: bcsc.c:389
pastix_parsec_init
void pastix_parsec_init(pastix_data_t *pastix, int *argc, char **argv[], const int *bindtab)
Startup the PaRSEC runtime system.
Definition: parsec.c:59
solverBackupRestore
int solverBackupRestore(SolverMatrix *solvmtx, const SolverBackup_t *b)
Restore initial values.
Definition: solver_backup.c:140
DPARM_EPSILON_MAGN_CTRL
@ DPARM_EPSILON_MAGN_CTRL
Definition: api.h:156
PastixSchedStarPU
@ PastixSchedStarPU
Definition: api.h:314
pastix_clrcores.h
kernelsTraceStop
double kernelsTraceStop(const pastix_data_t *pastix_data)
Stop the trace module.
Definition: kernels_trace.c:160
PASTIX_ERR_BADPARAMETER
@ PASTIX_ERR_BADPARAMETER
Definition: api.h:351
pastix_task_numfact
int pastix_task_numfact(pastix_data_t *pastix_data, spmatrix_t *spm)
Perform all the numerical factorization steps: fill the internal block CSC and the solver matrix stru...
Definition: pastix_task_sopalin.c:587
PastixFactLDLH
@ PastixFactLDLH
Definition: api.h:304
IPARM_FACTORIZATION
@ IPARM_FACTORIZATION
Definition: api.h:99
pastix_subtask_sopalin
int pastix_subtask_sopalin(pastix_data_t *pastix_data)
Factorize the given problem using Cholesky or LU decomposition.
Definition: pastix_task_sopalin.c:376
IPARM_COMPRESS_PRESELECT
@ IPARM_COMPRESS_PRESELECT
Definition: api.h:132