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  spmatrix_t *spmg = NULL;
89  double time;
90 
91  /*
92  * Check parameters
93  */
94  if (pastix_data == NULL) {
95  pastix_print_error( "pastix_subtask_spm2bcsc: wrong pastix_data parameter" );
97  }
98  if (spm == NULL) {
99  pastix_print_error( "pastix_subtask_spm2bcsc: wrong spm parameter" );
101  }
102  if ( !(pastix_data->steps & STEP_ANALYSE) ) {
103  pastix_print_error( "pastix_subtask_spm2bcsc: All steps from pastix_task_init() to pastix_task_blend() have to be called before calling this function" );
105  }
106 
107  /*
108  * Compute the norm of A, to scale the epsilon parameter for pivoting
109  */
110  {
111  pastix_data->dparm[ DPARM_A_NORM ] = spmNorm( SpmFrobeniusNorm, spm );
112  if ( pastix_data->iparm[IPARM_VERBOSE] > PastixVerboseNo ) {
113  pastix_print( pastix_data->inter_node_procnum, 0,
114  " ||A||_2 = %e\n",
115  pastix_data->dparm[ DPARM_A_NORM ] );
116  }
117  }
118 
119  /*
120  * Switch the SolverMatrix pointer according to the runtime.
121  */
122  {
123  pastix_int_t isched = pastix_data->iparm[IPARM_SCHEDULER];
124 
125  /*
126  * Temporary disable LDL[th] factorization with runtime
127  */
128  if ( ( isSchedRuntime( isched ) ) &&
129  ( ( pastix_data->iparm[IPARM_FACTORIZATION] == PastixFactLDLH ) ||
130  ( pastix_data->iparm[IPARM_FACTORIZATION] == PastixFactLDLT ) ) ) {
131  isched = PastixSchedDynamic;
132  }
133 
134  if ( isSchedRuntime( isched ) ) {
135  pastix_data->solvmatr = pastix_data->solvglob;
136  }
137  else {
138  assert( isSchedPthread( isched ) );
139  pastix_data->solvmatr = pastix_data->solvloc;
140  }
141  pastix_data->sched = isched;
142  }
143 
144  /*
145  * Fill in the internal blocked CSC. We consider that if this step is called
146  * the spm values have changed so we need to update the blocked csc.
147  */
148  if ( pastix_data->bcsc != NULL )
149  {
150  bcscExit( pastix_data->bcsc );
151  memFree_null( pastix_data->bcsc );
152  }
153 
154  MALLOC_INTERN( pastix_data->bcsc, 1, pastix_bcsc_t );
155 
156  if ( spm->loc2glob == NULL ) {
157  spmg = spm;
158  }
159 #if defined(PASTIX_WITH_MPI)
160  else {
161  if ( pastix_data->iparm[IPARM_VERBOSE] > PastixVerboseNot ) {
162  pastix_print(pastix_data->procnum, 0, "bcscInit: the SPM has to be centralized for the moment\n");
163  }
164  spmg = malloc( sizeof(spmatrix_t) );
165  spmGather( spm, -1, spmg );
166  }
167 #endif
168 
169  time = bcscInit( spmg,
170  pastix_data->ordemesh,
171  pastix_data->solvmatr,
172  (pastix_data->iparm[IPARM_FACTORIZATION] == PastixFactLU), /*&& (! pastix_data->iparm[IPARM_ONLY_REFINE]) )*/
173  pastix_data->bcsc );
174 
175  /* Free the gathered spm */
176  if ( spmg != spm ) {
177  spmExit( spmg );
178  memFree_null( spmg );
179  }
180 
181  if ( pastix_data->iparm[IPARM_VERBOSE] > PastixVerboseNot ) {
182  pastix_print( pastix_data->inter_node_procnum, 0, OUT_BCSC_TIME, time );
183  }
184 
185  if ( pastix_data->iparm[IPARM_FREE_CSCUSER] ) {
186  spmExit( spm );
187  }
188 
189  /*
190  * Invalidate following step, and add current step to the ones performed
191  */
192  pastix_data->steps &= ~STEP_BCSC2CTAB;
193  pastix_data->steps |= STEP_CSC2BCSC;
194 
195  return PASTIX_SUCCESS;
196 }
197 
198 /**
199  *******************************************************************************
200  *
201  * @ingroup pastix_numfact
202  *
203  * @brief Fill the internal solver matrix structure.
204  *
205  * This step is linked with the pastix_subtask_sopalin() since this structure is
206  * only used during the numerical factorization.
207  *
208  * This routine is affected by the following parameters:
209  * IPARM_VERBOSE, [IPARM_FACTORIZATION.
210  *
211  *******************************************************************************
212  *
213  * @param[inout] pastix_data
214  * The pastix_data structure that describes the solver instance.
215  * On exit, the internal solver structure is filled with entries from
216  * the internal block CSC matrix.
217  *
218  *******************************************************************************
219  *
220  * @retval PASTIX_SUCCESS on successful exit,
221  * @retval PASTIX_ERR_BADPARAMETER if one parameter is incorrect,
222  * @retval PASTIX_ERR_OUTOFMEMORY if one allocation failed.
223  *
224  *******************************************************************************/
225 int
226 pastix_subtask_bcsc2ctab( pastix_data_t *pastix_data )
227 {
228  pastix_bcsc_t *bcsc;
229  pastix_lr_t *lr;
230  Clock timer;
231  int mtxtype;
232 
233  /*
234  * Check parameters
235  */
236  if (pastix_data == NULL) {
237  pastix_print_error( "pastix_subtask_bcsc2ctab: wrong pastix_data parameter" );
239  }
240  if ( !(pastix_data->steps & STEP_CSC2BCSC) ) {
241  pastix_print_error( "pastix_subtask_bcsc2ctab: All steps from pastix_task_init() to pastix_stask_blend() have to be called before calling this function" );
243  }
244  if (pastix_data->bcsc == NULL) {
245  pastix_print_error( "pastix_subtask_bcsc2ctab: wrong pastix_data->bcsc parameter" );
247  }
248 
249  /* Ensure that the scheduler is correct */
250  pastix_check_and_correct_scheduler( pastix_data );
251 
252  clockSyncStart( timer, pastix_data->inter_node_comm );
253 
254  /* Initialize low-rank parameters */
255  lr = &(pastix_data->solvmatr->lowrank);
256  lr->compress_when = pastix_data->iparm[IPARM_COMPRESS_WHEN];
257  lr->compress_method = pastix_data->iparm[IPARM_COMPRESS_METHOD];
258  lr->compress_min_width = pastix_data->iparm[IPARM_COMPRESS_MIN_WIDTH];
259  lr->compress_min_height = pastix_data->iparm[IPARM_COMPRESS_MIN_HEIGHT];
260  lr->compress_preselect = pastix_data->iparm[IPARM_COMPRESS_PRESELECT];
261  lr->use_reltol = pastix_data->iparm[IPARM_COMPRESS_RELTOL];
262  lr->tolerance = pastix_data->dparm[DPARM_COMPRESS_TOLERANCE];
263  lr->ilu_lvl = pastix_data->iparm[IPARM_COMPRESS_ILUK];
264 
265  pastix_lr_minratio = pastix_data->dparm[DPARM_COMPRESS_MIN_RATIO];
266  pastix_lr_ortho = pastix_data->iparm[IPARM_COMPRESS_ORTHO];
267 
268  bcsc = pastix_data->bcsc;
269  lr->core_ge2lr = ge2lrMethods[ pastix_data->iparm[IPARM_COMPRESS_METHOD] ][bcsc->flttype-2];
270  lr->core_rradd = rraddMethods[ pastix_data->iparm[IPARM_COMPRESS_METHOD] ][bcsc->flttype-2];
271 
272  /*
273  * Update the ilu level based on when
274  */
275  if ( lr->ilu_lvl == -2 ) {
276  if ( pastix_data->iparm[IPARM_COMPRESS_WHEN] == PastixCompressWhenBegin ) {
277  lr->ilu_lvl = -1;
278  }
279  else if ( pastix_data->iparm[IPARM_COMPRESS_WHEN] == PastixCompressWhenEnd ) {
280  lr->ilu_lvl = INT_MAX;
281  }
282  else if ( pastix_data->iparm[IPARM_COMPRESS_WHEN] == PastixCompressWhenDuring ) {
283  lr->ilu_lvl = 0;
284  }
285  }
286 
287  /*
288  * Set the maximum rank function
289  */
290  if ( pastix_data->iparm[IPARM_COMPRESS_WHEN] == PastixCompressWhenBegin ) {
292  }
293  else {
295  }
296 
297  pastix_data->solvmatr->factotype = pastix_data->iparm[IPARM_FACTORIZATION];
298 
299 #if defined(PASTIX_WITH_MPI)
300  if ( ( pastix_data->iparm[IPARM_COMPRESS_WHEN] != PastixCompressNever ) &&
301  ( pastix_data->iparm[IPARM_SCHEDULER] == PastixSchedParsec ) &&
302  ( pastix_data->procnbr > 1 ) )
303  {
304  pastix_print_error( "pastix_task_sopalin: Low-Rank with MPI communication is not available yet with PaRSEC\n" );
306  }
307 #endif
308  /*
309  * Fill in the internal coeftab structure. We consider that if this step is
310  * called the bcsc values have changed, or a factorization have already been
311  * performed, so we need to update the coeftab arrays.
312  */
313  if (pastix_data->bcsc != NULL)
314  {
315  coeftabExit( pastix_data->solvmatr );
316  }
317 
318  coeftabInit( pastix_data,
319  pastix_data->iparm[IPARM_FACTORIZATION] == PastixFactLU ? PastixLUCoef : PastixLCoef );
320 
321  switch( pastix_data->iparm[IPARM_FACTORIZATION] ) {
322  case PastixFactLLH:
323  case PastixFactLDLH:
324  mtxtype = PastixHermitian;
325  break;
326 
327  case PastixFactLLT:
328  case PastixFactLDLT:
329  mtxtype = PastixHermitian;
330  break;
331 
332  case PastixFactLU:
333  default:
334  mtxtype = PastixGeneral;
335  }
336 
337 #if defined(PASTIX_WITH_PARSEC)
338  if ( pastix_data->iparm[IPARM_SCHEDULER] == PastixSchedParsec )
339  {
340  /* Start PaRSEC if not already started */
341  if (pastix_data->parsec == NULL) {
342  int argc = 0;
343  pastix_parsec_init( pastix_data, &argc, NULL, NULL );
344  }
345  /* Create the matrix descriptor */
346  parsec_sparse_matrix_init( pastix_data->solvmatr,
347  pastix_size_of( bcsc->flttype ), mtxtype,
348  pastix_data->inter_node_procnbr,
349  pastix_data->inter_node_procnum );
350  }
351 #endif
352 
353 #if defined(PASTIX_WITH_STARPU)
354  if ( pastix_data->iparm[IPARM_SCHEDULER] == PastixSchedStarPU )
355  {
356  /* Create the matrix descriptor */
357  starpu_sparse_matrix_init( pastix_data->solvmatr,
358  mtxtype,
359  pastix_data->inter_node_procnbr,
360  pastix_data->inter_node_procnum,
361  pastix_data->bcsc->flttype );
362  }
363 #endif
364 
365  clockSyncStop( timer, pastix_data->inter_node_comm );
366  if ( pastix_data->iparm[IPARM_VERBOSE] > PastixVerboseNot ) {
367  pastix_print( pastix_data->inter_node_procnum, 0, OUT_COEFTAB_TIME,
368  clockVal(timer) );
369  }
370 
371  /* Invalidate following step, and add current step to the ones performed */
372  pastix_data->steps &= ~STEP_NUMFACT;
373  pastix_data->steps |= STEP_BCSC2CTAB;
374 
375  (void)mtxtype;
376  return PASTIX_SUCCESS;
377 }
378 
379 /**
380  *******************************************************************************
381  *
382  * @ingroup pastix_numfact
383  *
384  * @brief Factorize the given problem using Cholesky or LU decomposition.
385  *
386  * The user can call the pastix_task_solve() to obtain the solution.
387  *
388  * This routine is affected by the following parameters:
389  * IPARM_VERBOSE, IPARM_FACTORIZATION, IPARM_COMPRESS_WHEN.
390  *
391  *******************************************************************************
392  *
393  * @param[inout] pastix_data
394  * The pastix_data structure that describes the solver instance.
395  * On exit, the solver matrix structure stores the factorization of the
396  * given problem.
397  *
398  *******************************************************************************
399  *
400  * @retval PASTIX_SUCCESS on successful exit,
401  * @retval PASTIX_ERR_BADPARAMETER if one parameter is incorrect,
402  * @retval PASTIX_ERR_OUTOFMEMORY if one allocation failed.
403  *
404  *******************************************************************************/
405 int
406 pastix_subtask_sopalin( pastix_data_t *pastix_data )
407 {
408  sopalin_data_t sopalin_data;
409  SolverBackup_t *sbackup;
410  pastix_bcsc_t *bcsc;
411  PASTIX_Comm pastix_comm;
412  pastix_int_t *iparm;
413  double *dparm;
414 
415  /*
416  * Check parameters
417  */
418  if (pastix_data == NULL) {
419  pastix_print_error( "pastix_subtask_sopalin: wrong pastix_data parameter" );
421  }
422  if ( !(pastix_data->steps & STEP_ANALYSE) ) {
423  pastix_print_error( "pastix_subtask_sopalin: All steps from pastix_task_init() to pastix_task_analyze() have to be called before calling this function" );
425  }
426  if ( !(pastix_data->steps & STEP_CSC2BCSC) ) {
427  pastix_print_error( "pastix_subtask_sopalin: All steps from pastix_task_init() to pastix_task_analyze() have to be called before calling this function" );
429  }
430  if ( !(pastix_data->steps & STEP_BCSC2CTAB) ) {
431  pastix_print_error( "pastix_subtask_sopalin: All steps from pastix_task_init() to pastix_task_analyze() have to be called before calling this function" );
433  }
434 
435  bcsc = pastix_data->bcsc;
436  if (bcsc == NULL) {
437  pastix_print_error( "pastix_subtask_sopalin: wrong pastix_data_bcsc parameter" );
439  }
440 
441  /* Ensure that the scheduler is correct */
442  pastix_check_and_correct_scheduler( pastix_data );
443 
444  pastix_comm = pastix_data->inter_node_comm;
445  iparm = pastix_data->iparm;
446  dparm = pastix_data->dparm;
447 
448  /* Prepare the sopalin_data structure */
449  {
450  double threshold;
451 
452  sopalin_data.solvmtx = pastix_data->solvmatr;
453 
454  /* TODO: might change the behavior: if the user wants a ratio of the norm, it could compute it himself */
455  if ( dparm[ DPARM_EPSILON_MAGN_CTRL ] < 0. ) {
456  threshold = - dparm[ DPARM_EPSILON_MAGN_CTRL ];
457  }
458  else if ( dparm[ DPARM_EPSILON_MAGN_CTRL ] == 0. ) {
459  /*
460  * Use the rule presented in "Making Sparse Gaussian Elimination
461  * Scalable by Static Pivoting", S. X. Li, J. Demmel
462  *
463  * sqrt(eps) * ||A||_1 is a small half precision perturbation
464  * introduced in the pb when the pivot is too small
465  */
466  double eps;
467  if ( (bcsc->flttype == PastixFloat) || (bcsc->flttype == PastixComplex32) ) {
468  eps = LAPACKE_slamch_work( 'e' );
469  }
470  else {
471  eps = LAPACKE_dlamch_work( 'e' );
472  }
473  threshold = sqrt(eps) * dparm[DPARM_A_NORM];
474  }
475  else {
476  threshold = dparm[ DPARM_EPSILON_MAGN_CTRL ] * dparm[DPARM_A_NORM];
477  }
478 
479  sopalin_data.solvmtx->diagthreshold = threshold;
480  sopalin_data.solvmtx->nbpivots = 0;
481 
482  sopalin_data.cpu_coefs = &(pastix_data->cpu_models->coefficients[bcsc->flttype-2]);
483  sopalin_data.gpu_coefs = &(pastix_data->gpu_models->coefficients[bcsc->flttype-2]);
484  sopalin_data.cpu_models = pastix_data->cpu_models;
485  sopalin_data.gpu_models = pastix_data->gpu_models;
486  }
487 
488  sbackup = solverBackupInit( pastix_data->solvmatr );
489  pastix_data->solvmatr->restore = 2;
490  {
491  void (*factofct)( pastix_data_t *, sopalin_data_t *);
492  double timer, timer_local, flops, flops_g = 0.;
493 
494  factofct = sopalinFacto[ pastix_data->iparm[IPARM_FACTORIZATION] ][bcsc->flttype-2];
495  assert(factofct);
496 
497  kernelsTraceStart( pastix_data );
498  clockSyncStart( timer, pastix_comm );
499  clockStart(timer_local);
500 
501  factofct( pastix_data, &sopalin_data );
502 
503  clockStop(timer_local);
504  clockSyncStop( timer, pastix_comm );
505  kernelsTraceStop( pastix_data );
506 
507  /* Output time and flops */
508  pastix_data->dparm[DPARM_FACT_TIME] = clockVal(timer);
509  flops = pastix_data->dparm[DPARM_FACT_THFLOPS] / pastix_data->dparm[DPARM_FACT_TIME];
510  pastix_data->dparm[DPARM_FACT_FLOPS] = ((flops / 1024.) / 1024.) / 1024.;
511 
512  /* Compute the number of flops performed (may be different from THFLOPS for low-rank) */
513  {
514  double flops_l = pastix_data->dparm[DPARM_FACT_RLFLOPS];
515  MPI_Reduce( &flops_l, &flops_g, 1, MPI_DOUBLE, MPI_SUM, 0, pastix_data->inter_node_comm );
516  }
517 
518  /* Reduce the number of pivots */
519  {
520  int nbpivot_l = sopalin_data.solvmtx->nbpivots;
521  int nbpivot_g;
522  MPI_Allreduce( &nbpivot_l, &nbpivot_g, 1, MPI_INT, MPI_SUM, pastix_data->inter_node_comm );
523  pastix_data->iparm[IPARM_STATIC_PIVOTING] = nbpivot_g;
524  }
525 
526  if ( iparm[IPARM_VERBOSE] > PastixVerboseNot ) {
527  pastix_print( pastix_data->inter_node_procnum, 0, OUT_SOPALIN_TIME,
528  clockVal(timer),
529  pastix_print_value( flops ),
530  pastix_print_unit( flops ),
531  pastix_print_value( flops_g ),
532  pastix_print_unit( flops_g ),
533  (long)pastix_data->iparm[IPARM_STATIC_PIVOTING] );
534  }
535 
536 #if defined(PASTIX_WITH_PARSEC) && defined(PASTIX_DEBUG_PARSEC)
537  {
538  int i;
539  fprintf(stderr, "-- Check status of PaRSEC nbtasks counters\n" );
540  for (i=0; i<32; i++) {
541  if ( parsec_nbtasks_on_gpu[i] != 0 ) {
542  fprintf(stderr, "Device %d => %d\n",
543  i, parsec_nbtasks_on_gpu[i] );
544  }
545  parsec_nbtasks_on_gpu[i] = 0;
546  }
547  }
548 #endif /* defined(PASTIX_WITH_PARSEC) && defined(PASTIX_DEBUG_PARSEC) */
549  }
550  solverBackupRestore( pastix_data->solvmatr, sbackup );
551  solverBackupExit( sbackup );
552 
553 #if defined(PASTIX_NUMFACT_DUMP_SOLVER)
554  {
555  FILE *stream = NULL;
556 
557  pastix_gendirectories( pastix_data );
558  stream = pastix_fopenw( pastix_data->dir_local, "solver.eps", "w" );
559  if ( stream ) {
560  solverDraw( pastix_data->solvmatr,
561  stream, iparm[IPARM_VERBOSE],
562  pastix_data->dir_local );
563  fclose(stream);
564  }
565  }
566 #endif
567 
568  if ( (pastix_data->iparm[IPARM_VERBOSE] > PastixVerboseNot) &&
569  (pastix_data->iparm[IPARM_COMPRESS_WHEN] != PastixCompressNever) )
570  {
571  /* Compute the memory gain */
572  coeftabMemory[bcsc->flttype-2]( pastix_data->solvmatr, pastix_data->dparm );
573 #if defined(PASTIX_WITH_MPI)
574  MPI_Allreduce( MPI_IN_PLACE, pastix_data->dparm + DPARM_MEM_FR, 1, MPI_DOUBLE, MPI_SUM, pastix_data->inter_node_comm );
575  MPI_Allreduce( MPI_IN_PLACE, pastix_data->dparm + DPARM_MEM_LR, 1, MPI_DOUBLE, MPI_SUM, pastix_data->inter_node_comm );
576 #endif
577  }
578 
579  /* Invalidate following steps, and add factorization step to the ones performed */
580  pastix_data->steps &= ~( STEP_BCSC2CTAB |
581  STEP_SOLVE |
582  STEP_REFINE );
583  pastix_data->steps |= STEP_NUMFACT;
584 
585  (void)pastix_comm;
586  return EXIT_SUCCESS;
587 }
588 
589 /**
590  *******************************************************************************
591  *
592  * @ingroup pastix_users
593  *
594  * @brief Perform all the numerical factorization steps:
595  * fill the internal block CSC and the solver matrix structures,
596  * then apply the factorization step.
597  *
598  * The user can call the pastix_task_solve() to obtain the solution.
599  *
600  * This routine is affected by the following parameters:
601  * IPARM_VERBOSE.
602  *
603  *******************************************************************************
604  *
605  * @param[inout] pastix_data
606  * The pastix_data structure that describes the solver instance.
607  * On exit, the internal block CSC is filled with entries from the
608  * spm matrix and the solver matrix structure stores the factorization
609  * of the given problem.
610  *
611  * @param[inout] spm
612  * The sparse matrix descriptor that describes problem instance.
613  * On exit, the spm structure is freed if IPARM_FREE_CSCUSER is true.
614  *
615  *******************************************************************************
616  *
617  * @retval PASTIX_SUCCESS on successful exit,
618  * @retval PASTIX_ERR_BADPARAMETER if one parameter is incorrect,
619  * @retval PASTIX_ERR_OUTOFMEMORY if one allocation failed.
620  *
621  *******************************************************************************/
622 int
623 pastix_task_numfact( pastix_data_t *pastix_data,
624  spmatrix_t *spm )
625 {
626  pastix_int_t *iparm;
627  pastix_int_t procnum;
628  pastix_int_t steps;
629  int rc;
630 
631  /*
632  * Check parameters
633  */
634  if (pastix_data == NULL) {
635  pastix_print_error( "pastix_task_sopalin: wrong pastix_data parameter" );
637  }
638  if (spm == NULL) {
639  pastix_print_error( "pastix_task_sopalin: wrong spm parameter" );
641  }
642  if ( !(pastix_data->steps & STEP_ANALYSE) ) {
643  pastix_print_error( "pastix_task_sopalin: All steps from pastix_task_init() to pastix_task_blend() have to be called before calling this function" );
645  }
646 
647  iparm = pastix_data->iparm;
648  procnum = pastix_data->inter_node_procnum;
649 
650  if (iparm[IPARM_VERBOSE] > PastixVerboseNot) {
651  pastix_print( procnum, 0, OUT_STEP_SOPALIN,
652  pastixFactotypeStr( iparm[IPARM_FACTORIZATION] ) );
653  }
654 
655  /* Invalidate upcoming steps */
656  steps = ~(STEP_CSC2BCSC |
657  STEP_BCSC2CTAB |
658  STEP_NUMFACT |
659  STEP_SOLVE |
660  STEP_REFINE );
661  pastix_data->steps &= steps;
662 
663  rc = pastix_subtask_spm2bcsc( pastix_data, spm );
664  if (rc != PASTIX_SUCCESS) {
665  return rc;
666  }
667 
668  rc = pastix_subtask_bcsc2ctab( pastix_data );
669  if (rc != PASTIX_SUCCESS) {
670  return rc;
671  }
672 
673  rc = pastix_subtask_sopalin( pastix_data );
674  if (rc != PASTIX_SUCCESS) {
675  return rc;
676  }
677 
678  return EXIT_SUCCESS;
679 }
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:305
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:175
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:414
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
DPARM_MEM_FR
@ DPARM_MEM_FR
Definition: api.h:167
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
starpu_sparse_matrix_init
void starpu_sparse_matrix_init(SolverMatrix *solvmtx, int mtxtype, int nodes, int myrank, pastix_coeftype_t flttype)
Generate the StarPU descriptor of the sparse matrix.
Definition: starpu_sparse_matrix.c:159
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:364
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:435
DPARM_A_NORM
@ DPARM_A_NORM
Definition: api.h:174
IPARM_COMPRESS_MIN_HEIGHT
@ IPARM_COMPRESS_MIN_HEIGHT
Definition: api.h:127
DPARM_MEM_LR
@ DPARM_MEM_LR
Definition: api.h:168
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:315
pastix_lr_s::core_ge2lr
fct_ge2lr_t core_ge2lr
Definition: pastix_lowrank.h:157
PastixCompressWhenBegin
@ PastixCompressWhenBegin
Definition: api.h:365
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:302
PastixVerboseNot
@ PastixVerboseNot
Definition: api.h:209
PASTIX_SUCCESS
@ PASTIX_SUCCESS
Definition: api.h:346
PastixVerboseNo
@ PastixVerboseNo
Definition: api.h:210
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:456
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:605
PastixHermitian
@ PastixHermitian
Definition: api.h:437
DPARM_FACT_TIME
@ DPARM_FACT_TIME
Definition: api.h:163
DPARM_COMPRESS_MIN_RATIO
@ DPARM_COMPRESS_MIN_RATIO
Definition: api.h:176
IPARM_COMPRESS_ILUK
@ IPARM_COMPRESS_ILUK
Definition: api.h:133
coeftabExit
void coeftabExit(SolverMatrix *solvmtx)
Free the solver matrix structure.
Definition: coeftab.c:185
PastixCompressWhenDuring
@ PastixCompressWhenDuring
Definition: api.h:367
PastixFactLDLT
@ PastixFactLDLT
Definition: api.h:303
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:42
PastixLUCoef
@ PastixLUCoef
Definition: api.h:458
coeftab.h
PastixCompressWhenEnd
@ PastixCompressWhenEnd
Definition: api.h:366
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:304
pastix_subtask_bcsc2ctab
int pastix_subtask_bcsc2ctab(pastix_data_t *pastix_data)
Fill the internal solver matrix structure.
Definition: pastix_task_sopalin.c:226
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:138
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:384
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:316
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:353
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:623
PastixFactLDLH
@ PastixFactLDLH
Definition: api.h:306
IPARM_FACTORIZATION
@ IPARM_FACTORIZATION
Definition: api.h:99
PastixSchedDynamic
@ PastixSchedDynamic
Definition: api.h:317
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:406
IPARM_COMPRESS_PRESELECT
@ IPARM_COMPRESS_PRESELECT
Definition: api.h:132