PaStiX Handbook  6.2.1
pastix_task_solve.c
Go to the documentation of this file.
1 /**
2  *
3  * @file pastix_task_solve.c
4  *
5  * PaStiX solve 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 Theophile Terraz
16  * @author Tony Delarue
17  * @date 2021-07-02
18  *
19  **/
20 #define _GNU_SOURCE 1
21 #include "common.h"
22 #include "bcsc/bcsc.h"
23 #include "pastix/order.h"
24 #include "blend/solver.h"
25 #include "sopalin/sopalin_data.h"
26 
27 #include "bcsc/bcsc_z.h"
28 #include "bcsc/bcsc_c.h"
29 #include "bcsc/bcsc_d.h"
30 #include "bcsc/bcsc_s.h"
31 
32 #if defined(PASTIX_DEBUG_SOLVE)
33 #include <z_spm.h>
34 #include <c_spm.h>
35 #include <d_spm.h>
36 #include <s_spm.h>
37 
38 static volatile int32_t step = 0;
39 
40 static inline void
41 dump_rhs( pastix_data_t *pastix_data, const char *name,
42  spm_coeftype_t flttype, int m, int n, const void *b, int ldb )
43 {
44  FILE *f;
45  f = pastix_fopenw( (pastix_data->dir_local == NULL) ? "." : pastix_data->dir_local,
46  name, "w" );
47 
48  switch( flttype ) {
49  case SpmComplex64:
50  z_spmDensePrint( f, m, n, b, ldb );
51  break;
52  case SpmComplex32:
53  c_spmDensePrint( f, m, n, b, ldb );
54  break;
55  case SpmDouble:
56  d_spmDensePrint( f, m, n, b, ldb );
57  break;
58  case SpmFloat:
59  s_spmDensePrint( f, m, n, b, ldb );
60  break;
61  case SpmPattern:
62  ;
63  }
64 
65  fclose(f);
66 }
67 #else
68 static inline void
69 dump_rhs( pastix_data_t *pastix_data, const char *name,
70  spm_coeftype_t flttype, int m, int n, const void *b, int ldb )
71 {
72  (void)pastix_data;
73  (void)name;
74  (void)m;
75  (void)n;
76  (void)b;
77  (void)flttype;
78  (void)ldb;
79 }
80 #endif
81 
82 /**
83  *******************************************************************************
84  *
85  * @ingroup pastix_solve
86  *
87  * @brief Apply a permutation on the right-and-side vector before the solve step.
88  *
89  * This routine is affected by the following parameters:
90  * IPARM_VERBOSE, IPARM_FACTORIZATION, IPARM_APPLYPERM_WS.
91  *
92  *******************************************************************************
93  *
94  * @param[inout] pastix_data
95  * The pastix_data structure that describes the solver instance.
96  *
97  * @param[in] flttype
98  * This arithmetic of the sparse matrix.
99  *
100  * @param[in] dir
101  * Forward or backword application of the permutation.
102  *
103  * @param[in] m
104  * Size of the right-and-side vectors.
105  *
106  * @param[in] n
107  * Number of right-and-side vectors.
108  *
109  * @param[inout] b
110  * The right-and-side vectors (can be multiple RHS).
111  *
112  * @param[in] ldb
113  * The leading dimension of the right-and-side vectors.
114  *
115  *******************************************************************************
116  *
117  * @retval PASTIX_SUCCESS on successful exit,
118  * @retval PASTIX_ERR_BADPARAMETER if one parameter is incorrect.
119  *
120  *******************************************************************************/
121 int
122 pastix_subtask_applyorder( pastix_data_t *pastix_data,
123  pastix_coeftype_t flttype, pastix_dir_t dir,
124  pastix_int_t m, pastix_int_t n, void *b, pastix_int_t ldb )
125 {
126  pastix_int_t *perm;
127  int ts;
128 
129  /*
130  * Check parameters
131  */
132  if (pastix_data == NULL) {
133  errorPrint("pastix_subtask_applyorder: wrong pastix_data parameter");
135  }
136  if (b == NULL) {
137  errorPrint("pastix_subtask_applyorder: wrong b parameter");
139  }
140  if ( !(pastix_data->steps & STEP_CSC2BCSC) ) {
141  errorPrint("pastix_subtask_applyorder: All steps from pastix_task_init() to pastix_subtask_csc2bcsc() have to be called before calling this function");
143  }
144 
145  /* Make sure ordering is 0 based */
146  if ( pastix_data->ordemesh->baseval != 0 ) {
147  errorPrint("pastix_subtask_applyorder: ordermesh must be 0-based");
149  }
150 
151  ts = pastix_data->iparm[IPARM_APPLYPERM_WS];
152  perm = pastix_data->ordemesh->peritab;
153 
154  /* See also xlapmr and xlapmt */
155  switch( flttype ) {
156  case PastixComplex64:
157  bvec_zlapmr( ts, dir, m, n, b, ldb, perm );
158  break;
159 
160  case PastixComplex32:
161  bvec_clapmr( ts, dir, m, n, b, ldb, perm );
162  break;
163 
164  case PastixFloat:
165  bvec_slapmr( ts, dir, m, n, b, ldb, perm );
166  break;
167 
168  case PastixDouble:
169  default:
170  bvec_dlapmr( ts, dir, m, n, b, ldb, perm );
171  }
172 
173 #if defined(PASTIX_DEBUG_SOLVE)
174  {
175  char *fullname;
176  int32_t lstep = pastix_atomic_inc_32b( &step );
177  asprintf( &fullname, "solve_%02d_applyorder_%s.rhs", lstep,
178  ( dir == PastixDirForward ) ? "Forward" : "Backward" );
179 
180  dump_rhs( pastix_data, fullname,
181  flttype, m, n, b, ldb );
182  free(fullname);
183  }
184 #endif
185 
186  return PASTIX_SUCCESS;
187 }
188 
189 /**
190  *******************************************************************************
191  *
192  * @ingroup pastix_solve
193  *
194  * @brief Apply a triangular solve on the right-and-side vectors.
195  *
196  * This routine is affected by the following parameters:
197  * IPARM_VERBOSE, IPARM_FACTORIZATION.
198  *
199  *******************************************************************************
200  *
201  * @param[inout] pastix_data
202  * The pastix_data structure that describes the solver instance.
203  *
204  * @param[in] flttype
205  * This arithmetic of the sparse matrix.
206  *
207  * @param[in] side
208  * Left or right application.
209  *
210  * @param[in] uplo
211  * Upper or Lower part.
212  *
213  * @param[in] trans
214  * With or without transposition (or conjugate transposition).
215  *
216  * @param[in] diag
217  * Diagonal terms are unit or not.
218  *
219  * @param[in] nrhs
220  * The number of right-and-side vectors.
221  *
222  * @param[inout] b
223  * The right-and-side vector (can be multiple RHS).
224  * On exit, the solution is stored in place of the right-hand-side vector.
225  *
226  * @param[in] ldb
227  * The leading dimension of the right-and-side vectors.
228  *
229  *******************************************************************************
230  *
231  * @retval PASTIX_SUCCESS on successful exit,
232  * @retval PASTIX_ERR_BADPARAMETER if one parameter is incorrect.
233  *
234  *******************************************************************************/
235 int
236 pastix_subtask_trsm( pastix_data_t *pastix_data,
237  pastix_coeftype_t flttype, pastix_side_t side,
238  pastix_uplo_t uplo, pastix_trans_t trans, pastix_diag_t diag,
239  pastix_int_t nrhs, void *b, pastix_int_t ldb )
240 {
241  sopalin_data_t sopalin_data;
242  int i, bs = nrhs;
243 
244 #if defined(PASTIX_WITH_MPI)
245  bs = 1;
246 #endif
247 
248  /*
249  * Check parameters
250  */
251  if (pastix_data == NULL) {
252  errorPrint("pastix_subtask_trsm: wrong pastix_data parameter");
254  }
255  if (b == NULL) {
256  errorPrint("pastix_subtask_trsm: wrong b parameter");
258  }
259  if ( !(pastix_data->steps & STEP_NUMFACT) ) {
260  errorPrint("pastix_subtask_trsm: All steps from pastix_task_init() to pastix_task_numfact() have to be called before calling this function");
262  }
263 
264  /*
265  * Ensure that the scheduler is correct and is in the same
266  * family that the one used for the factorization.
267  */
268  pastix_check_and_correct_scheduler( pastix_data );
269 
270  sopalin_data.solvmtx = pastix_data->solvmatr;
271 
272  switch (flttype) {
273  case PastixComplex64:
274  {
275  pastix_complex64_t *lb = b;
276  for( i = 0; i < nrhs; i+=bs, lb += ldb ) {
277  sopalin_ztrsm( pastix_data, side, uplo, trans, diag,
278  &sopalin_data, bs, lb, ldb );
279  }
280  }
281  break;
282  case PastixComplex32:
283  {
284  pastix_complex32_t *lb = b;
285  for( i = 0; i < nrhs; i+=bs, lb += ldb ) {
286  sopalin_ctrsm( pastix_data, side, uplo, trans, diag,
287  &sopalin_data, bs, lb, ldb );
288  }
289  }
290  break;
291  case PastixDouble:
292  {
293  double *lb = b;
294  trans = (trans == PastixConjTrans) ? PastixTrans : trans;
295  for( i = 0; i < nrhs; i+=bs, lb += ldb ) {
296  sopalin_dtrsm( pastix_data, side, uplo, trans, diag,
297  &sopalin_data, bs, lb, ldb );
298  }
299  }
300  break;
301  case PastixFloat:
302  {
303  float *lb = b;
304  trans = (trans == PastixConjTrans) ? PastixTrans : trans;
305  for( i = 0; i < nrhs; i+=bs, lb += ldb ) {
306  sopalin_strsm( pastix_data, side, uplo, trans, diag,
307  &sopalin_data, bs, lb, ldb );
308  }
309  }
310  break;
311  default:
312  fprintf(stderr, "Unknown floating point arithmetic\n" );
313  }
314 
315 #if defined(PASTIX_DEBUG_SOLVE)
316  {
317  char *fullname;
318  int32_t lstep = pastix_atomic_inc_32b( &step );
319  asprintf( &fullname, "solve_%02d_trsm_%c%c%c%c.rhs", lstep,
320  (side == PastixLeft) ? 'L' : 'R',
321  (uplo == PastixUpper) ? 'U' : 'L',
322  (trans == PastixConjTrans) ? 'C' : (trans == PastixTrans ? 'T' : 'N'),
323  (diag == PastixUnit) ? 'U' : 'N' );
324 
325  dump_rhs( pastix_data, fullname,
326  flttype, pastix_data->bcsc->gN, nrhs, b, ldb );
327  free(fullname);
328  }
329 #endif
330 
331  return PASTIX_SUCCESS;
332 }
333 
334 /**
335  *******************************************************************************
336  *
337  * @ingroup pastix_solve
338  *
339  * @brief Apply a diagonal operation on the right-and-side vectors.
340  *
341  * This routine is affected by the following parameters:
342  * IPARM_VERBOSE, IPARM_FACTORIZATION.
343  *
344  *******************************************************************************
345  *
346  * @param[inout] pastix_data
347  * The pastix_data structure that describes the solver instance.
348  *
349  * @param[in] flttype
350  * This arithmetic of the sparse matrix.
351  *
352  * @param[in] nrhs
353  * The number of right-and-side vectors.
354  *
355  * @param[inout] b
356  * The right-and-side vector (can be multiple RHS).
357  * On exit, the solution is stored in place of the right-hand-side vector.
358  *
359  * @param[in] ldb
360  * The leading dimension of the right-and-side vectors.
361  *
362  *******************************************************************************
363  *
364  * @retval PASTIX_SUCCESS on successful exit,
365  * @retval PASTIX_ERR_BADPARAMETER if one parameter is incorrect.
366  *
367  *******************************************************************************/
368 int
369 pastix_subtask_diag( pastix_data_t *pastix_data, pastix_coeftype_t flttype,
370  pastix_int_t nrhs, void *b, pastix_int_t ldb )
371 {
372  sopalin_data_t sopalin_data;
373 
374  /*
375  * Check parameters
376  */
377  if (pastix_data == NULL) {
378  errorPrint("pastix_subtask_diag: wrong pastix_data parameter");
380  }
381  if (b == NULL) {
382  errorPrint("pastix_subtask_diag: wrong b parameter");
384  }
385  if ( !(pastix_data->steps & STEP_NUMFACT) ) {
386  errorPrint("pastix_subtask_trsm: All steps from pastix_task_init() to pastix_task_numfact() have to be called before calling this function");
388  }
389 
390  /*
391  * Ensure that the scheduler is correct and is in the same
392  * family that the one used for the factorization.
393  */
394  pastix_check_and_correct_scheduler( pastix_data );
395 
396  sopalin_data.solvmtx = pastix_data->solvmatr;
397 
398  switch (flttype) {
399  case PastixComplex64:
400  sopalin_zdiag( pastix_data, &sopalin_data, nrhs, (pastix_complex64_t *)b, ldb );
401  break;
402  case PastixComplex32:
403  sopalin_cdiag( pastix_data, &sopalin_data, nrhs, (pastix_complex32_t *)b, ldb );
404  break;
405  case PastixDouble:
406  sopalin_ddiag( pastix_data, &sopalin_data, nrhs, (double *)b, ldb );
407  break;
408  case PastixFloat:
409  sopalin_sdiag( pastix_data, &sopalin_data, nrhs, (float *)b, ldb );
410  break;
411  default:
412  fprintf(stderr, "Unknown floating point arithmetic\n" );
413  }
414 
415 #if defined(PASTIX_DEBUG_SOLVE)
416  {
417  char *fullname;
418  int32_t lstep = pastix_atomic_inc_32b( &step );
419  asprintf( &fullname, "solve_%02d_diag.rhs", lstep );
420 
421  dump_rhs( pastix_data, fullname,
422  flttype, pastix_data->bcsc->gN, nrhs, b, ldb );
423  free(fullname);
424  }
425 #endif
426 
427  return PASTIX_SUCCESS;
428 }
429 
430 /**
431  *******************************************************************************
432  *
433  * @ingroup pastix_solve
434  *
435  * @brief Solve the given problem without applying the permutation.
436  *
437  * @warning The input vector is considered already permuted. For a solve step
438  * with permutation, see pastix_task_solve()
439  *
440  * This routine is affected by the following parameters:
441  * IPARM_VERBOSE, IPARM_FACTORIZATION.
442  *
443  *******************************************************************************
444  *
445  * @param[inout] pastix_data
446  * The pastix_data structure that describes the solver instance.
447  *
448  * @param[in] transA
449  * PastixNoTrans: A is not transposed (CSC matrix)
450  * PastixTrans: A is transposed (CSR of symmetric/general matrix)
451  * PastixConjTrans: A is conjugate transposed (CSR of hermitian matrix)
452  *
453  * @param[in] nrhs
454  * The number of right-and-side vectors.
455  *
456  * @param[inout] b
457  * The right-and-side vectors (can be multiple RHS).
458  * On exit, the solution is stored in place of the right-hand-side vector.
459  *
460  * @param[in] ldb
461  * The leading dimension of the right-and-side vectors.
462  *
463  *******************************************************************************
464  *
465  * @retval PASTIX_SUCCESS on successful exit,
466  * @retval PASTIX_ERR_BADPARAMETER if one parameter is incorrect.
467  *
468  *******************************************************************************/
469 int
470 pastix_subtask_solve_adv( pastix_data_t *pastix_data, pastix_trans_t transA,
471  pastix_int_t nrhs, void *b, pastix_int_t ldb )
472 {
473  pastix_bcsc_t *bcsc;
474  pastix_factotype_t factotype;
475 
476  pastix_uplo_t uplo;
477  pastix_diag_t diag;
478  pastix_trans_t trans;
479  pastix_trans_t transfact = PastixTrans;
480 
481  /*
482  * Check parameters
483  */
484  if (pastix_data == NULL) {
485  errorPrint("pastix_task_solve: wrong pastix_data parameter");
487  }
488  if ( !(pastix_data->steps & STEP_NUMFACT) ) {
489  errorPrint("pastix_task_solve: All steps from pastix_task_init() to pastix_task_numfact() have to be called before calling this function");
491  }
492 
493  bcsc = pastix_data->bcsc;
494  factotype = pastix_data->iparm[IPARM_FACTORIZATION];
495 
496  /**
497  *
498  * Summary of the operations to perform the full solve if the original A was
499  * either transposed or conjugate transposed.
500  *
501  * op(A) | Factorization | Step 1 | Step 2
502  * ----------+---------------------+-----------+-----------
503  * NoTrans | L U | L y = b | U x = y
504  * NoTrans | L L^t | L y = b | L^t x = y
505  * NoTrans | L L^h | L y = b | L^h x = y
506  * Trans |(L U )^t = U^t L^t | U^t y = b | L^t x = y
507  * Trans |(L L^t)^t = L L^t | L y = b | L^t x = y
508  * Trans |(L L^h)^t = c(L) L^t | Not handled (c(L))
509  * ConjTrans |(L U )^h = U^h L^h | Not handled (U^h)
510  * ConjTrans |(L L^t)^h = c(L) L^h | Not handled (c(L))
511  * ConjTrans |(L L^h)^h = L L^h | L y = b | L^h x = y
512  *
513  */
514  /* Value of the transpose case */
515  if ( ((bcsc->flttype == PastixComplex32) || (bcsc->flttype == PastixComplex64)) &&
516  ((factotype == PastixFactLLH ) || ( factotype == PastixFactLDLH )) )
517  {
518  transfact = PastixConjTrans;
519  }
520 
521  if ( (transA != PastixNoTrans) &&
522  (transA != transfact) )
523  {
524  errorPrint("pastix_task_solve: transA incompatible with the factotype used (require extra conj(L) not handled)");
526  }
527 
528  {
529  double timer;
530  /* Start timer */
531  clockSyncStart( timer, pastix_data->inter_node_comm );
532 
533  /*
534  * Solve the first step
535  */
536  uplo = PastixLower;
537  trans = PastixNoTrans;
538  diag = PastixNonUnit;
539 
540  if (( transA != PastixNoTrans ) && ( factotype == PastixFactLU )) {
541  uplo = PastixUpper;
542  trans = transA;
543  }
544  if (( transA == PastixNoTrans ) && ( factotype == PastixFactLU )) {
545  diag = PastixUnit;
546  }
547 
548  if (( factotype == PastixFactLDLT ) ||
549  ( factotype == PastixFactLDLH ) )
550  {
551  diag = PastixUnit;
552  }
553 
554  pastix_subtask_trsm( pastix_data, bcsc->flttype,
555  PastixLeft, uplo, trans, diag,
556  nrhs, b, ldb );
557 
558  /*
559  * Solve the diagonal step
560  */
561  if( (factotype == PastixFactLDLT) ||
562  (factotype == PastixFactLDLH) )
563  {
564  /* Solve y = D z with z = ([L^t | L^h] P x) */
565  pastix_subtask_diag( pastix_data, bcsc->flttype, nrhs, b, ldb );
566  }
567 
568  /*
569  * Solve the second step
570  */
571  uplo = PastixLower;
572  trans = transfact;
573  diag = PastixNonUnit;
574 
575  if (( transA == PastixNoTrans ) && ( factotype == PastixFactLU )) {
576  uplo = PastixUpper;
577  trans = PastixNoTrans;
578  }
579  if (( transA != PastixNoTrans ) && ( factotype == PastixFactLU )) {
580  diag = PastixUnit;
581  }
582 
583  if (( factotype == PastixFactLDLT ) ||
584  ( factotype == PastixFactLDLH ) )
585  {
586  diag = PastixUnit;
587  }
588 
589  pastix_subtask_trsm( pastix_data, bcsc->flttype,
590  PastixLeft, uplo, trans, diag,
591  nrhs, b, ldb );
592 
593  /* Stop Timer */
594  clockSyncStop( timer, pastix_data->inter_node_comm );
595 
596  pastix_data->dparm[DPARM_SOLV_TIME] = clockVal(timer);
597  if ( pastix_data->iparm[IPARM_VERBOSE] > PastixVerboseNot ) {
598  pastix_print( pastix_data->inter_node_procnum, 0, OUT_TIME_SOLV,
599  pastix_data->dparm[DPARM_SOLV_TIME] );
600  }
601  }
602 
603  return EXIT_SUCCESS;
604 }
605 
606 /**
607  *******************************************************************************
608  *
609  * @ingroup pastix_solve
610  *
611  * @brief Solve the given problem without applying the permutation.
612  *
613  * @warning The input vector is considered already permuted. For a solve step
614  * with permutation, see pastix_task_solve()
615  *
616  * This routine is affected by the following parameters:
617  * IPARM_VERBOSE, IPARM_FACTORIZATION, IPARM_TRANSPOSE_SOLVE.
618  *
619  *******************************************************************************
620  *
621  * @param[inout] pastix_data
622  * The pastix_data structure that describes the solver instance.
623  *
624  * @param[in] nrhs
625  * The number of right-and-side vectors.
626  *
627  * @param[inout] b
628  * The right-and-side vectors (can be multiple RHS).
629  * On exit, the solution is stored in place of the right-hand-side vector.
630  *
631  * @param[in] ldb
632  * The leading dimension of the right-and-side vectors.
633  *
634  *******************************************************************************
635  *
636  * @retval PASTIX_SUCCESS on successful exit,
637  * @retval PASTIX_ERR_BADPARAMETER if one parameter is incorrect.
638  *
639  *******************************************************************************/
640 int
641 pastix_subtask_solve( pastix_data_t *pastix_data,
642  pastix_int_t nrhs, void *b, pastix_int_t ldb )
643 {
644  return pastix_subtask_solve_adv( pastix_data,
645  pastix_data->iparm[IPARM_TRANSPOSE_SOLVE],
646  nrhs, b, ldb );
647 }
648 
649 /**
650  *******************************************************************************
651  *
652  * @ingroup pastix_users
653  *
654  * @brief Solve the given problem.
655  *
656  * This routine is affected by the following parameters:
657  * IPARM_VERBOSE, IPARM_FACTORIZATION, IPARM_TRANSPOSE_SOLVE.
658  *
659  *******************************************************************************
660  *
661  * @param[inout] pastix_data
662  * The pastix_data structure that describes the solver instance.
663  *
664  * @param[in] nrhs
665  * The number of right-and-side vectors.
666  *
667  * @param[inout] b
668  * The right-and-side vectors (can be multiple RHS).
669  * On exit, the solution is stored in place of the right-hand-side vector.
670  *
671  * @param[in] ldb
672  * The leading dimension of the right-and-side vectors.
673  *
674  *******************************************************************************
675  *
676  * @retval PASTIX_SUCCESS on successful exit,
677  * @retval PASTIX_ERR_BADPARAMETER if one parameter is incorrect.
678  *
679  *******************************************************************************/
680 int
681 pastix_task_solve( pastix_data_t *pastix_data,
682  pastix_int_t nrhs, void *b, pastix_int_t ldb )
683 {
684  pastix_bcsc_t *bcsc;
685  void *bglob = NULL;
686  void *tmp = NULL;
687 
688  /*
689  * Check parameters
690  */
691  if (pastix_data == NULL) {
692  errorPrint("pastix_task_solve: wrong pastix_data parameter");
694  }
695 
696  if ( !(pastix_data->steps & STEP_NUMFACT) ) {
697  errorPrint("pastix_task_solve: Numerical factorization hasn't been done.");
699  }
700 
701  bcsc = pastix_data->bcsc;
702 
703  /* The spm is distributed, so we have to gather the RHS */
704  if ( pastix_data->csc->loc2glob != NULL ) {
705  if( pastix_data->iparm[IPARM_VERBOSE] > PastixVerboseNo ) {
706  pastix_print( pastix_data->procnum, 0, "pastix_task_solve: the RHS has to be centralized for the moment\n" );
707  }
708 
709  tmp = b;
710  spmGatherRHS( nrhs, pastix_data->csc, b, ldb, &bglob, -1 );
711  b = bglob;
712  ldb = pastix_data->csc->gNexp;
713  }
714 
715  /* Compute P * b */
716  pastix_subtask_applyorder( pastix_data, bcsc->flttype,
717  PastixDirForward, bcsc->gN, nrhs, b, ldb );
718 
719  /* Solve A x = b */
720  pastix_subtask_solve( pastix_data, nrhs, b, ldb );
721 
722  /* Compute P^t * b */
723  pastix_subtask_applyorder( pastix_data, bcsc->flttype,
724  PastixDirBackward, bcsc->gN, nrhs, b, ldb );
725 
726  if( tmp != NULL ) {
727  pastix_int_t ldbglob = ldb;
728  ldb = pastix_data->csc->nexp;
729  b = tmp;
730 
731  spmExtractLocalRHS( nrhs, pastix_data->csc, bglob, ldbglob, b, ldb );
732  memFree_null( bglob );
733  }
734 
735  return EXIT_SUCCESS;
736 }
solver.h
PastixTrans
@ PastixTrans
Definition: api.h:423
pastix_uplo_t
enum pastix_uplo_e pastix_uplo_t
Upper/Lower part.
pastix_subtask_applyorder
int pastix_subtask_applyorder(pastix_data_t *pastix_data, pastix_coeftype_t flttype, pastix_dir_t dir, pastix_int_t m, pastix_int_t n, void *b, pastix_int_t ldb)
Apply a permutation on the right-and-side vector before the solve step.
Definition: pastix_task_solve.c:122
bcsc_s.h
pastix_bcsc_s
Internal column block distributed CSC matrix.
Definition: bcsc.h:37
bvec_zlapmr
int bvec_zlapmr(int thread_safe, pastix_dir_t dir, pastix_int_t m, pastix_int_t n, pastix_complex64_t *A, pastix_int_t lda, pastix_int_t *perm)
Apply a row permutation to a matrix A (LAPACK xlatmr)
Definition: bvec_zcompute.c:878
step
bvec_slapmr
int bvec_slapmr(int thread_safe, pastix_dir_t dir, pastix_int_t m, pastix_int_t n, float *A, pastix_int_t lda, pastix_int_t *perm)
Apply a row permutation to a matrix A (LAPACK xlatmr)
Definition: bvec_scompute.c:878
bvec_dlapmr
int bvec_dlapmr(int thread_safe, pastix_dir_t dir, pastix_int_t m, pastix_int_t n, double *A, pastix_int_t lda, pastix_int_t *perm)
Apply a row permutation to a matrix A (LAPACK xlatmr)
Definition: bvec_dcompute.c:878
pastix_trans_t
enum pastix_trans_e pastix_trans_t
Transpostion.
pastix_subtask_diag
int pastix_subtask_diag(pastix_data_t *pastix_data, pastix_coeftype_t flttype, pastix_int_t nrhs, void *b, pastix_int_t ldb)
Apply a diagonal operation on the right-and-side vectors.
Definition: pastix_task_solve.c:369
PastixLower
@ PastixLower
Definition: api.h:443
pastix_task_solve
int pastix_task_solve(pastix_data_t *pastix_data, pastix_int_t nrhs, void *b, pastix_int_t ldb)
Solve the given problem.
Definition: pastix_task_solve.c:681
PastixDirForward
@ PastixDirForward
Definition: api.h:489
PastixNoTrans
@ PastixNoTrans
Definition: api.h:422
PastixUpper
@ PastixUpper
Definition: api.h:442
PastixConjTrans
@ PastixConjTrans
Definition: api.h:424
pastix_bcsc_s::flttype
int flttype
Definition: bcsc.h:41
PastixNonUnit
@ PastixNonUnit
Definition: api.h:463
bcsc.h
pastix_factotype_t
enum pastix_factotype_e pastix_factotype_t
Factorization algorithms available for IPARM_FACTORIZATION parameter.
PastixFactLLH
@ PastixFactLLH
Definition: api.h:300
IPARM_TRANSPOSE_SOLVE
@ IPARM_TRANSPOSE_SOLVE
Definition: api.h:105
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
PastixDirBackward
@ PastixDirBackward
Definition: api.h:490
IPARM_APPLYPERM_WS
@ IPARM_APPLYPERM_WS
Definition: api.h:107
bcsc_d.h
PastixUnit
@ PastixUnit
Definition: api.h:464
pastix_dir_t
enum pastix_dir_e pastix_dir_t
Direction.
bcsc_c.h
bvec_clapmr
int bvec_clapmr(int thread_safe, pastix_dir_t dir, pastix_int_t m, pastix_int_t n, pastix_complex32_t *A, pastix_int_t lda, pastix_int_t *perm)
Apply a row permutation to a matrix A (LAPACK xlatmr)
Definition: bvec_ccompute.c:878
PastixFactLDLT
@ PastixFactLDLT
Definition: api.h:301
pastix_diag_t
enum pastix_diag_e pastix_diag_t
Diagonal.
DPARM_SOLV_TIME
@ DPARM_SOLV_TIME
Definition: api.h:167
IPARM_VERBOSE
@ IPARM_VERBOSE
Definition: api.h:36
PastixFactLU
@ PastixFactLU
Definition: api.h:302
order.h
bcsc_z.h
pastix_coeftype_t
#define pastix_coeftype_t
Arithmetic types.
Definition: api.h:281
pastix_subtask_trsm
int pastix_subtask_trsm(pastix_data_t *pastix_data, pastix_coeftype_t flttype, pastix_side_t side, pastix_uplo_t uplo, pastix_trans_t trans, pastix_diag_t diag, pastix_int_t nrhs, void *b, pastix_int_t ldb)
Apply a triangular solve on the right-and-side vectors.
Definition: pastix_task_solve.c:236
pastix_bcsc_s::gN
int gN
Definition: bcsc.h:38
pastix_subtask_solve_adv
int pastix_subtask_solve_adv(pastix_data_t *pastix_data, pastix_trans_t transA, pastix_int_t nrhs, void *b, pastix_int_t ldb)
Solve the given problem without applying the permutation.
Definition: pastix_task_solve.c:470
PastixLeft
@ PastixLeft
Definition: api.h:471
PASTIX_ERR_BADPARAMETER
@ PASTIX_ERR_BADPARAMETER
Definition: api.h:351
pastix_side_t
enum pastix_side_e pastix_side_t
Side of the operation.
PastixFactLDLH
@ PastixFactLDLH
Definition: api.h:304
IPARM_FACTORIZATION
@ IPARM_FACTORIZATION
Definition: api.h:99
pastix_subtask_solve
int pastix_subtask_solve(pastix_data_t *pastix_data, pastix_int_t nrhs, void *b, pastix_int_t ldb)
Solve the given problem without applying the permutation.
Definition: pastix_task_solve.c:641