PaStiX Handbook  6.3.2
parsec_sparse_matrix.c
Go to the documentation of this file.
1 /**
2  *
3  * @file parsec_sparse_matrix.c
4  *
5  * PaStiX sparse matrix descriptor for parsec.
6  *
7  * @copyright 2016-2023 Bordeaux INP, CNRS (LaBRI UMR 5800), Inria,
8  * Univ. Bordeaux. All rights reserved.
9  *
10  * @version 6.3.2
11  * @author Mathieu Faverge
12  * @author Pierre Ramet
13  * @date 2023-07-21
14  *
15  * @addtogroup pastix_parsec
16  * @{
17  *
18  **/
19 #include "common.h"
20 #include <parsec.h>
21 #include <parsec/data.h>
22 #include <parsec/data_distribution.h>
23 #if defined(PASTIX_CUDA_FERMI)
24 #include <parsec/devices/cuda/dev_cuda.h>
25 #endif
26 
27 #include "blend/solver.h"
28 #include "pastix_parsec.h"
29 
30 /**
31  *******************************************************************************
32  *
33  * @brief Compute the triplet (uplo, cblknum, bloknum) from the key.
34  *
35  * This function convert the unique key identifier of each piece of data, to its
36  * original information.
37  *
38  *******************************************************************************
39  *
40  * @param[in] key
41  * The key identifier to decode.
42  *
43  * @param[in] solvmtx
44  * The solver matrix structure to access information about the blocks
45  * and cblks.
46  *
47  * @param[out] uplo
48  * On exit, uplo is 0 if the key corresponds to the lower part, or 1
49  * for the upper part.
50  *
51  * @param[out] cblknum
52  * On exit, the index of the cblk encoded in the key.
53  *
54  * @param[out] bloknum
55  * On exit, the index of the blok encoded in the key.
56  * bloknum = 0, if the piece of data is the full cblk for 1D kernels.
57  * bloknum > 0, if the piece of data is the (bloknum-1)^th block in the cblk.
58  *
59  ******************************************************************************/
60 static inline void
61 spm_data_key_to_value( parsec_data_key_t key,
62  const SolverMatrix *solvmtx,
63  int *uplo,
64  pastix_int_t *cblknum,
65  pastix_int_t *bloknum )
66 {
67  parsec_data_key_t key2;
68  pastix_int_t cblkmin2d, cblknbr;
69 
70  cblknbr = solvmtx->cblknbr;
71  cblkmin2d = solvmtx->cblkmin2d;
72  key2 = 2 * cblknbr;
73 
74  /* This is a block */
75  if ( key >= key2 ) {
76  pastix_int_t m, n, ld;
77 
78  key2 = key - key2;
79  ld = solvmtx->cblkmaxblk * 2;
80 
81  m = key2 % ld;
82  n = key2 / ld;
83 
84  *uplo = m % 2;
85  *bloknum = m / 2;
86  *cblknum = cblkmin2d + n;
87  }
88  /* This is a cblk */
89  else {
90  *uplo = key % 2;
91  *cblknum = key / 2;
92  *bloknum = -1;
93  }
94 }
95 
96 /**
97  *******************************************************************************
98  *
99  * @brief Compute the unique key from the triplet (uplo, cblknum, bloknum).
100  *
101  * This function convert the unique triplet in a unique key identifier for each
102  * piece of data.
103  *
104  *******************************************************************************
105  *
106  * @param[in] mat
107  * The sparse matrix descriptor.
108  *
109  * @param[in] ...
110  * This function must receive three int:
111  * - uplo: 0 for the lower part, 1 for the upper part
112  * - cblknum: the index of the cblk
113  * - bloknum: =0 if this is the full cblk,
114  * >0 for the (bloknum-1)^th block in the cblk.
115  *
116  *******************************************************************************
117  *
118  * @return The unique key identifier for the piece of data.
119  *
120  ******************************************************************************/
121 static uint32_t
122 parsec_sparse_matrix_data_key( parsec_data_collection_t *mat, ... )
123 {
124  va_list ap;
126  int uplo;
127  pastix_int_t cblknum, bloknum;
128 
129  va_start(ap, mat);
130  uplo = va_arg(ap, int);
131  cblknum = va_arg(ap, int);
132  bloknum = va_arg(ap, int) - 1;
133  va_end(ap);
134 
135  uplo = uplo ? 1 : 0;
136 
137  if ( bloknum == -1 ) {
138  return cblknum * 2 + uplo;
139  }
140  else {
141  pastix_int_t offset, ld, cblknbr;
142  pastix_int_t cblkmin2d, n;
143 
144  cblknbr = spmtx->solvmtx->cblknbr;
145  cblkmin2d = spmtx->solvmtx->cblkmin2d;
146  ld = spmtx->solvmtx->cblkmaxblk * 2;
147  offset = cblknbr * 2;
148  n = cblknum - cblkmin2d;
149 
150  return offset + n * ld + bloknum * 2 + uplo;
151  }
152 }
153 
154 /**
155  *******************************************************************************
156  *
157  * @brief Return the rank of the owner of the piece of data (uplo, cblknum,
158  * bloknum).
159  *
160  *******************************************************************************
161  *
162  * @param[in] mat
163  * The sparse matrix descriptor.
164  *
165  * @param[in] ...
166  * This function must receive three int:
167  * - uplo: 0 for the lower part, 1 for the upper part
168  * - cblknum: the index of the cblk
169  * - bloknum: =0 if this is the full cblk,
170  * >0 for the (bloknum-1)^th block in the cblk.
171  *
172  *******************************************************************************
173  *
174  * @return The rank index of the owner of the data.
175  *
176  ******************************************************************************/
177 static uint32_t
178 parsec_sparse_matrix_rank_of( parsec_data_collection_t *mat, ... )
179 {
181  SolverCblk *cblk;
182  va_list ap;
183  int uplo;
184  pastix_int_t cblknum, bloknum;
185 
186  va_start(ap, mat);
187  uplo = va_arg(ap, int);
188  cblknum = va_arg(ap, int);
189  bloknum = va_arg(ap, int) - 1;
190  va_end(ap);
191 
192  uplo = uplo ? 1 : 0;
193 
194  cblk = spmtx->solvmtx->cblktab + cblknum;
195 
196  (void)bloknum;
197  return cblk->ownerid;
198 }
199 
200 /**
201  *******************************************************************************
202  *
203  * @brief Return the rank of the owner of the piece of data (key)
204  *
205  *******************************************************************************
206  *
207  * @param[in] mat
208  * The sparse matrix descriptor.
209  *
210  * @param[in] key
211  * The unique key idenifier of a piece of data.
212  *
213  *******************************************************************************
214  *
215  * @return The rank index of the owner of the data.
216  *
217  ******************************************************************************/
218 static uint32_t
219 parsec_sparse_matrix_rank_of_key( parsec_data_collection_t *mat,
220  parsec_data_key_t key )
221 {
223  SolverMatrix *solvmtx = spmtx->solvmtx;
224  SolverCblk *cblk;
225  int uplo;
226  pastix_int_t cblknum, bloknum;
227 
228  spm_data_key_to_value( key, solvmtx,
229  &uplo, &cblknum, &bloknum );
230 
231  cblk = solvmtx->cblktab + cblknum;
232 
233  return cblk->ownerid;
234 }
235 
236 /**
237  *******************************************************************************
238  *
239  * @brief Return the rank of the virtual process owner of the piece of data
240  * (uplo, cblknum, bloknum).
241  *
242  *******************************************************************************
243  *
244  * @param[in] mat
245  * The sparse matrix descriptor.
246  *
247  * @param[in] ...
248  * This function must receive three int:
249  * - uplo: 0 for the lower part, 1 for the upper part
250  * - cblknum: the index of the cblk
251  * - bloknum: =0 if this is the full cblk,
252  * >0 for the (bloknum-1)^th block in the cblk.
253  *
254  *******************************************************************************
255  *
256  * @return The rank index of the virtual process owner of the data.
257  *
258  ******************************************************************************/
259 static int32_t
260 parsec_sparse_matrix_vpid_of( parsec_data_collection_t *mat, ... )
261 {
262  (void)mat;
263  return 0;
264 }
265 
266 /**
267  *******************************************************************************
268  *
269  * @brief Return the rank of the virtual process owner of the piece of data (key)
270  *
271  *******************************************************************************
272  *
273  * @param[in] mat
274  * The sparse matrix descriptor.
275  *
276  * @param[in] key
277  * The unique key idenifier of a piece of data.
278  *
279  *******************************************************************************
280  *
281  * @return The rank index of the virtual process owner of the data.
282  *
283  ******************************************************************************/
284 static int32_t
285 parsec_sparse_matrix_vpid_of_key( parsec_data_collection_t *mat,
286  parsec_data_key_t key )
287 {
288  (void)mat; (void)key;
289  return 0;
290 }
291 
292 /**
293  *******************************************************************************
294  *
295  * @brief Return the data handler associated to the piece of data (uplo,
296  * cblknum, bloknum).
297  *
298  *******************************************************************************
299  *
300  * @param[in] mat
301  * The sparse matrix descriptor.
302  *
303  * @param[in] ...
304  * This function must receive three int:
305  * - uplo: 0 for the lower part, 1 for the upper part
306  * - cblknum: the index of the cblk
307  * - bloknum: =0 if this is the full cblk,
308  * >0 for the (bloknum-1)^th block in the cblk.
309  *
310  *******************************************************************************
311  *
312  * @return The pointer to the data handler of the data.
313  *
314  ******************************************************************************/
315 static parsec_data_t *
316 parsec_sparse_matrix_data_of( parsec_data_collection_t *mat, ... )
317 {
319  SolverCblk *cblk;
320  va_list ap;
321  int uplo;
322  pastix_int_t cblknum, bloknum;
323 
324  va_start(ap, mat);
325  uplo = va_arg(ap, int);
326  cblknum = va_arg(ap, int);
327  bloknum = va_arg(ap, int) - 1;
328  va_end(ap);
329 
330  uplo = uplo ? 1 : 0;
331 
332  cblk = spmtx->solvmtx->cblktab + cblknum;
333 
334  /* This is a cblk */
335  if ( bloknum == -1 ) {
336  assert( cblk->handler[uplo] );
337  return (parsec_data_t*)(cblk->handler[uplo]);
338  }
339  /* This is a blok */
340  else {
341  SolverBlok *blok = cblk->fblokptr + bloknum;
342 
343  assert( blok->handler[uplo] );
344  return (parsec_data_t*)(blok->handler[uplo]);
345  }
346 }
347 
348 /**
349  *******************************************************************************
350  *
351  * @brief Return the data handler associated to the piece of data (key).
352  *
353  *******************************************************************************
354  *
355  * @param[in] mat
356  * The sparse matrix descriptor.
357  *
358  * @param[in] key
359  * The unique key idenifier of a piece of data.
360  *
361  *******************************************************************************
362  *
363  * @return The pointer to the data handler of the data.
364  *
365  ******************************************************************************/
366 static parsec_data_t *
367 parsec_sparse_matrix_data_of_key( parsec_data_collection_t *mat,
368  parsec_data_key_t key )
369 {
371  SolverMatrix *solvmtx = spmtx->solvmtx;
372  SolverCblk *cblk;
373  int uplo;
374  pastix_int_t cblknum, bloknum;
375 
376  spm_data_key_to_value( key, solvmtx,
377  &uplo, &cblknum, &bloknum );
378 
379  cblk = solvmtx->cblktab + cblknum;
380 
381  /* This is a cblk */
382  if ( bloknum == -1 ) {
383  assert( cblk->handler[uplo] );
384  return (parsec_data_t*)(cblk->handler[uplo]);
385  }
386  /* This is a blok */
387  else {
388  SolverBlok *blok = cblk->fblokptr + bloknum;
389 
390  assert( blok->handler[uplo] );
391  return (parsec_data_t*)(blok->handler[uplo]);
392  }
393 }
394 
395 #if defined(PARSEC_PROF_TRACE)
396 /**
397  *******************************************************************************
398  *
399  * @brief Convert the uinque key identifier to a human readable string.
400  *
401  *******************************************************************************
402  *
403  * @param[in] mat
404  * The sparse matrix descriptor.
405  *
406  * @param[in] key
407  * The unique key idenifier of a piece of data.
408  *
409  * @param[inout] buffer
410  * The char buffer of size buffer_size that will receive the string
411  * describing the piece of data.
412  *
413  * @param[in] buffer_size
414  * The size of the buffer buffer. buffer_size > 0.
415  *
416  *******************************************************************************
417  *
418  * @return The number of characters printed (excluding the null byte used to end
419  * output to strings) as returned by snprintf.
420  *
421  ******************************************************************************/
422 static int
423 parsec_sparse_matrix_key_to_string( parsec_data_collection_t *mat,
424  uint32_t key,
425  char *buffer, uint32_t buffer_size )
426 {
428  int uplo;
429  pastix_int_t cblknum, bloknum;
430  int res;
431 
432  spm_data_key_to_value( key, spmtx->solvmtx,
433  &uplo, &cblknum, &bloknum );
434 
435  res = snprintf(buffer, buffer_size, "(%d, %ld, %ld)",
436  uplo, (long int)cblknum, (long int)bloknum);
437  if (res < 0)
438  {
439  printf("error in key_to_string for tile (%d, %ld, %ld) key: %u\n",
440  uplo, (long int)cblknum, (long int)bloknum, key);
441  }
442  return res;
443 }
444 #endif
445 
446 #if defined(PASTIX_CUDA_FERMI)
447 /**
448  *******************************************************************************
449  *
450  * @brief TODO
451  *******************************************************************************
452  *
453  * @param[in] spmtx
454  * TODO
455  *
456  * @param[in] solvmtx
457  * TODO
458  *
459  ******************************************************************************/
460 void
461 parsec_sparse_matrix_init_fermi( parsec_sparse_matrix_desc_t *spmtx,
462  const SolverMatrix *solvmtx )
463 {
464  gpu_device_t* gpu_device;
465  SolverBlok *blok;
466  pastix_int_t i, b, bloknbr, ndevices;
467  size_t size;
468  int *tmp, *bloktab;
469 
470  ndevices = parsec_devices_enabled();
471  if ( ndevices <= 2 )
472  return;
473 
474  bloknbr = solvmtx->bloknbr;
475  size = 2 * bloknbr * sizeof(int);
476 
477  /**
478  * Initialize array on CPU
479  */
480  bloktab = (int*)malloc( size );
481  tmp = bloktab;
482  for (b=0, blok = solvmtx->bloktab;
483  b < bloknbr;
484  b++, blok++, tmp+=2)
485  {
486  tmp[0] = blok->frownum;
487  tmp[1] = blok->lrownum;
488  }
489 
490  ndevices -= 2;
491  spmtx->gpu_blocktab = calloc(ndevices, sizeof(void*));
492 
493  fprintf(stderr, "ndevices = %ld\n", ndevices );
494  for(i = 0; i < ndevices; i++) {
495  if( NULL == (gpu_device = (gpu_device_t*)parsec_devices_get(i+2)) ) continue;
496 
497  fprintf(stderr, "cuda index = %d\n", gpu_device->cuda_index );
498  cudaSetDevice( gpu_device->cuda_index );
499 
500  cudaMalloc( &(spmtx->gpu_blocktab[i]),
501  size );
502 
503  cudaMemcpy( spmtx->gpu_blocktab[i],
504  bloktab, size,
505  cudaMemcpyHostToDevice );
506  }
507  free(bloktab);
508 }
509 
510 /**
511  *******************************************************************************
512  *
513  * @brief TODO
514  *******************************************************************************
515  *
516  * @param[in] spmtx
517  * TODO
518  *
519  ******************************************************************************/
520 void
521 parsec_sparse_matrix_destroy_fermi( parsec_sparse_matrix_desc_t *spmtx )
522 {
523  gpu_device_t* gpu_device;
524  pastix_int_t i, ndevices;
525 
526  ndevices = parsec_devices_enabled();
527  if ( ndevices <= 2 )
528  return;
529 
530  ndevices -= 2;
531  for(i = 0; i < ndevices; i++) {
532  if( NULL == (gpu_device = (gpu_device_t*)parsec_devices_get(i+2)) ) continue;
533 
534  cudaSetDevice( gpu_device->cuda_index );
535  cudaFree( spmtx->gpu_blocktab[i] );
536  }
537 
538  free( spmtx->gpu_blocktab );
539 }
540 #endif /*defined(PASTIX_CUDA_FERMI)*/
541 
542 /**
543  *******************************************************************************
544  *
545  * @brief TODO
546  *******************************************************************************
547  *
548  * @param[in] o
549  * TODO
550  *
551  * @param[in] handler
552  * TODO
553  *
554  * @param[in] id
555  * TODO
556  *
557  * @param[in] cblk
558  * TODO
559  *
560  * @param[in] side
561  * TODO
562  *
563  ******************************************************************************/
564 static inline void
565 pastix_parsec_register_cblk_lr( parsec_data_collection_t *o,
566  parsec_data_t **handler,
567  pastix_int_t id,
568  const SolverCblk *cblk,
569  int side )
570 {
571  pastix_lrblock_t *dataptr = cblk->fblokptr->LRblock[side];
572  size_t size = ( cblk[1].fblokptr - cblk[0].fblokptr ) * sizeof( pastix_lrblock_t );
573 
574  parsec_data_create( handler, o, id, dataptr, size );
575 }
576 
577 /**
578  *******************************************************************************
579  *
580  * @brief TODO
581  *******************************************************************************
582  *
583  * @param[in] o
584  * TODO
585  *
586  * @param[in] handler
587  * TODO
588  *
589  * @param[in] id
590  * TODO
591  *
592  * @param[in] spmtx
593  * TODO
594  *
595  * @param[in] cblk
596  * TODO
597  *
598  * @param[in] side
599  * TODO
600  *
601  ******************************************************************************/
602 static inline void
603 pastix_parsec_register_cblk_fr( parsec_data_collection_t *o,
604  parsec_data_t **handler,
605  pastix_int_t id,
606  const parsec_sparse_matrix_desc_t *spmtx,
607  const SolverCblk *cblk,
608  int side )
609 {
610  void *dataptr = ( side == PastixUCoef ) ? cblk->ucoeftab : cblk->lcoeftab;
611  size_t size = (size_t)cblk->stride * (size_t)cblk_colnbr( cblk ) * (size_t)spmtx->typesze;
612 
613  parsec_data_create( handler, o, id, dataptr, size );
614 }
615 
616 /**
617  *******************************************************************************
618  *
619  * @brief TODO
620  *******************************************************************************
621  *
622  * @param[in] o
623  * TODO
624  *
625  * @param[in] cblknum
626  * TODO
627  *
628  * @param[in] spmtx
629  * TODO
630  *
631  * @param[in] cblk
632  * TODO
633  *
634  ******************************************************************************/
635 static inline void
636 pastix_parsec_register_cblk( parsec_data_collection_t *o,
637  pastix_int_t cblknum,
638  const parsec_sparse_matrix_desc_t *spmtx,
639  const SolverCblk *cblk )
640 {
641  parsec_data_t **handler = (parsec_data_t **)( cblk->handler );
642 
643  if ( cblk->cblktype & CBLK_COMPRESSED ) {
644  pastix_parsec_register_cblk_lr( o, handler, cblknum * 2, cblk, PastixLCoef );
645 
646  if ( spmtx->mtxtype == PastixGeneral ) {
647  pastix_parsec_register_cblk_lr( o, handler + 1, cblknum * 2 + 1, cblk, PastixUCoef );
648  }
649  }
650  else {
651  pastix_parsec_register_cblk_fr( o, handler, cblknum * 2, spmtx, cblk, PastixLCoef );
652  if ( spmtx->mtxtype == PastixGeneral ) {
654  o, handler + 1, cblknum * 2 + 1, spmtx, cblk, PastixUCoef );
655  }
656  }
657 
658  assert(cblk->handler[0]);
659  assert(cblk->handler[1] || ( spmtx->mtxtype != PastixGeneral ));
660 }
661 
662 /**
663  *******************************************************************************
664  *
665  * @brief Generate the PaRSEC descriptor of the sparse matrix.
666  *
667  * This function creates the PaRSEC descriptor that will provide the data
668  * mapping and memory location to PaRSEC for the computation.
669  *
670  *******************************************************************************
671  *
672  * @param[inout] solvmtx
673  * The solver matrix structure that describes the sparse matrix for
674  * PaStiX.
675  *
676  * @param[in] typesize
677  * The memory size of the arithmetic used to store the matrix
678  * coefficients.
679  *
680  * @param[in] mtxtype
681  * The type of sparse matrix to describe.
682  * @arg PastixGeneral: The sparse matrix is general.
683  * @arg PastixSymmetric: The sparse matrix is lower triangular symmetric.
684  * @arg PastixHermitian: The sparse matrix is lower triangular hermitian.
685  *
686  * @param[in] nodes
687  * The number of processes used to solve the problem.
688  *
689  * @param[in] myrank
690  * The rank of the calling process.
691  *
692  ******************************************************************************/
693 void
695  int typesize, pastix_mtxtype_t mtxtype,
696  int nodes, int myrank )
697 {
698  parsec_sparse_matrix_desc_t *spmtx = solvmtx->parsec_desc;
699  parsec_data_collection_t *o;
700  pastix_int_t cblknbr, cblkmin2d, ld;
701  parsec_data_key_t key1, key2;
702  SolverCblk *cblk;
703  SolverBlok *blok, *fblok, *lblok;
704  pastix_int_t m=0, n=0, cblknum, nbrow;
705  size_t size, offset;
706  char *ptrL, *ptrU;
707  pastix_lrblock_t *dataptrL, *dataptrU;
708 
709  if ( spmtx != NULL ) {
711  }
712  else {
714  }
715 
716  o = (parsec_data_collection_t*)spmtx;
717  parsec_data_collection_init( o, nodes, myrank );
718 
719  o->data_key = parsec_sparse_matrix_data_key;
720 #if defined(PARSEC_PROF_TRACE)
721  o->key_to_string = parsec_sparse_matrix_key_to_string;
722 #endif
723 
724  o->rank_of = parsec_sparse_matrix_rank_of;
725  o->rank_of_key = parsec_sparse_matrix_rank_of_key;
726  o->vpid_of = parsec_sparse_matrix_vpid_of;
727  o->vpid_of_key = parsec_sparse_matrix_vpid_of_key;
728  o->data_of = parsec_sparse_matrix_data_of;
729  o->data_of_key = parsec_sparse_matrix_data_of_key;
730 
731  spmtx->typesze = typesize;
732  spmtx->mtxtype = mtxtype;
733  spmtx->solvmtx = solvmtx;
734 
735  cblknbr = solvmtx->cblknbr;
736  cblkmin2d = solvmtx->cblkmin2d;
737  ld = solvmtx->cblkmaxblk * 2;
738  key1 = 2 * cblknbr;
739 
740  /* Initialize 1D cblk handlers */
741  cblk = spmtx->solvmtx->cblktab;
742  for(cblknum = 0;
743  cblknum < cblkmin2d;
744  cblknum++, n++, cblk++ )
745  {
746  if ( cblk->ownerid != myrank ) {
747  continue;
748  } else {
749  pastix_parsec_register_cblk( o, cblknum, spmtx, cblk );
750  }
751  }
752 
753  /* Initialize 2D cblk handlers */
754  cblk = spmtx->solvmtx->cblktab + cblkmin2d;
755  for(cblknum = cblkmin2d, n = 0;
756  cblknum < cblknbr;
757  cblknum++, n++, cblk++ )
758  {
759  if ( cblk->ownerid != myrank ) {
760  continue;
761  }
762 
763  pastix_parsec_register_cblk( o, cblknum, spmtx, cblk );
764 
765  if ( !(cblk->cblktype & CBLK_TASKS_2D) ) {
766  continue;
767  }
768 
769  if ( cblk->cblktype & CBLK_COMPRESSED ) {
770  /*
771  * Diagonal block
772  */
773  dataptrL = cblk->fblokptr->LRblock[0];
774  dataptrU = cblk->fblokptr->LRblock[1];
775  blok = cblk->fblokptr;
776  size = sizeof( pastix_lrblock_t );
777  offset = 0;
778  key2 = n * ld;
779 
780  assert( offset == 0 );
781  parsec_data_create( (parsec_data_t **)&( blok->handler[0] ),
782  o, key1 + key2, dataptrL + offset, size );
783 
784  if ( mtxtype == PastixGeneral ) {
785  parsec_data_create( (parsec_data_t **)&( blok->handler[1] ),
786  o, key1 + key2 + 1, dataptrU + offset, size );
787  }
788  else {
789  blok->handler[1] = NULL;
790  }
791 
792  /*
793  * Off-diagonal blocks
794  */
795  blok++;
796  key2 += 2;
797  lblok = cblk[1].fblokptr;
798  for ( ; blok < lblok; blok++, key2 += 2 ) {
799  fblok = blok;
800  m = 0;
801  size = blok_rownbr( blok );
802  offset = blok - cblk->fblokptr;
803  nbrow = 1;
804 
805  while ( ( blok < lblok ) && ( blok[0].fcblknm == blok[1].fcblknm ) &&
806  ( blok[0].lcblknm == blok[1].lcblknm ) ) {
807  blok++;
808  m++;
809  nbrow++;
810  }
811  size = nbrow;
812 
813  parsec_data_create( (parsec_data_t **)&( fblok->handler[0] ),
814  &spmtx->super,
815  key1 + key2,
816  dataptrL + offset,
817  size );
818 
819  if ( mtxtype == PastixGeneral ) {
820  parsec_data_create( (parsec_data_t **)&( fblok->handler[1] ),
821  &spmtx->super,
822  key1 + key2 + 1,
823  dataptrU + offset,
824  size );
825  }
826  else {
827  fblok->handler[1] = NULL;
828  }
829 
830  key2 += m * 2;
831  }
832  }
833  else {
834  /*
835  * Diagonal block
836  */
837  ptrL = cblk->lcoeftab;
838  ptrU = cblk->ucoeftab;
839  blok = cblk->fblokptr;
840  size = blok_rownbr( blok ) * cblk_colnbr( cblk ) * (size_t)spmtx->typesze;
841  offset = blok->coefind * (size_t)spmtx->typesze;
842  key2 = n * ld;
843 
844  assert( offset == 0 );
845  parsec_data_create(
846  (parsec_data_t **)&( blok->handler[0] ), o, key1 + key2, ptrL + offset, size );
847 
848  if ( mtxtype == PastixGeneral ) {
849  parsec_data_create( (parsec_data_t **)&( blok->handler[1] ),
850  o,
851  key1 + key2 + 1,
852  ptrU + offset,
853  size );
854  }
855  else {
856  blok->handler[1] = NULL;
857  }
858 
859  /*
860  * Off-diagonal blocks
861  */
862  blok++;
863  key2 += 2;
864  lblok = cblk[1].fblokptr;
865  for ( ; blok < lblok; blok++, key2 += 2 ) {
866  fblok = blok;
867  m = 0;
868  size = blok_rownbr( blok );
869  offset = blok->coefind * (size_t)spmtx->typesze;
870 
871  while ( ( blok < lblok ) &&
872  ( blok[0].fcblknm == blok[1].fcblknm ) &&
873  ( blok[0].lcblknm == blok[1].lcblknm ) )
874  {
875  blok++;
876  m++;
877  size += blok_rownbr( blok );
878  }
879  size *= cblk_colnbr( cblk ) * (size_t)spmtx->typesze;
880 
881  parsec_data_create( (parsec_data_t **)&( fblok->handler[0] ),
882  &spmtx->super,
883  key1 + key2,
884  ptrL + offset,
885  size );
886 
887  if ( mtxtype == PastixGeneral ) {
888  parsec_data_create( (parsec_data_t **)&( fblok->handler[1] ),
889  &spmtx->super,
890  key1 + key2 + 1,
891  ptrU + offset,
892  size );
893  }
894  else {
895  fblok->handler[1] = NULL;
896  }
897 
898  key2 += m * 2;
899  }
900  }
901  }
902 
903 #if defined( PASTIX_CUDA_FERMI )
904  parsec_sparse_matrix_init_fermi( spmtx, solvmtx );
905 #endif
906 
907  solvmtx->parsec_desc = spmtx;
908 }
909 
910 /**
911  *******************************************************************************
912  *
913  * @brief Free the PaRSEC descriptor of the sparse matrix.
914  *
915  * This function destroys the PaRSEC descriptor, but do not free the matrix data
916  * that are managed by PaStiX.
917  *
918  *******************************************************************************
919  *
920  * @param[inout] spmtx
921  * The descriptor to free.
922  *
923  ******************************************************************************/
924 void
926 {
927  SolverCblk *cblk;
928  SolverBlok *blok;
929  pastix_int_t i, cblkmin2d;
930 
931 #if defined(PASTIX_CUDA_FERMI)
932  parsec_sparse_matrix_destroy_fermi( spmtx );
933 #endif
934 
935  cblkmin2d = spmtx->solvmtx->cblkmin2d;
936  cblk = spmtx->solvmtx->cblktab;
937  for(i=0; i<cblkmin2d; i++, cblk++)
938  {
939  if ( cblk->handler[0] ) {
940  parsec_data_destroy( cblk->handler[0] );
941 
942  if ( spmtx->mtxtype == PastixGeneral ) {
943  parsec_data_destroy( cblk->handler[1] );
944  }
945  }
946 
947  cblk->handler[0] = NULL;
948  cblk->handler[1] = NULL;
949  }
950 
951  for(i=cblkmin2d; i<spmtx->solvmtx->cblknbr; i++, cblk++)
952  {
953  if ( cblk->handler[0] ) {
954  parsec_data_destroy( cblk->handler[0] );
955  if ( spmtx->mtxtype == PastixGeneral ) {
956  parsec_data_destroy( cblk->handler[1] );
957  }
958  }
959 
960  cblk->handler[0] = NULL;
961  cblk->handler[1] = NULL;
962 
963  blok = cblk->fblokptr;
964  while( blok < cblk[1].fblokptr )
965  {
966  if ( blok->handler[0] ) {
967  parsec_data_destroy( blok->handler[0] );
968  if ( spmtx->mtxtype == PastixGeneral ) {
969  parsec_data_destroy( blok->handler[1] );
970  }
971  }
972 
973  blok->handler[0] = NULL;
974  blok->handler[1] = NULL;
975 
976  blok++;
977  }
978  }
979 
980  parsec_data_collection_destroy( (parsec_data_collection_t*)spmtx );
981 }
982 
983 /**
984  *@}
985  */
BEGIN_C_DECLS typedef int pastix_int_t
Definition: datatypes.h:51
struct pastix_lrblock_s pastix_lrblock_t
The block low-rank structure to hold a matrix in low-rank form.
The block low-rank structure to hold a matrix in low-rank form.
#define PastixGeneral
Definition: api.h:458
spm_mtxtype_t pastix_mtxtype_t
Matrix symmetry type property.
Definition: api.h:457
@ PastixLCoef
Definition: api.h:478
@ PastixUCoef
Definition: api.h:479
parsec_data_collection_t super
Definition: pastix_parsec.h:36
pastix_mtxtype_t mtxtype
Definition: pastix_parsec.h:38
void parsec_sparse_matrix_init(SolverMatrix *solvmtx, int typesize, pastix_mtxtype_t mtxtype, int nodes, int myrank)
Generate the PaRSEC descriptor of the sparse matrix.
static uint32_t parsec_sparse_matrix_data_key(parsec_data_collection_t *mat,...)
Compute the unique key from the triplet (uplo, cblknum, bloknum).
static void pastix_parsec_register_cblk_lr(parsec_data_collection_t *o, parsec_data_t **handler, pastix_int_t id, const SolverCblk *cblk, int side)
TODO.
static int32_t parsec_sparse_matrix_vpid_of_key(parsec_data_collection_t *mat, parsec_data_key_t key)
Return the rank of the virtual process owner of the piece of data (key)
static void pastix_parsec_register_cblk(parsec_data_collection_t *o, pastix_int_t cblknum, const parsec_sparse_matrix_desc_t *spmtx, const SolverCblk *cblk)
TODO.
static parsec_data_t * parsec_sparse_matrix_data_of_key(parsec_data_collection_t *mat, parsec_data_key_t key)
Return the data handler associated to the piece of data (key).
static uint32_t parsec_sparse_matrix_rank_of(parsec_data_collection_t *mat,...)
Return the rank of the owner of the piece of data (uplo, cblknum, bloknum).
static void pastix_parsec_register_cblk_fr(parsec_data_collection_t *o, parsec_data_t **handler, pastix_int_t id, const parsec_sparse_matrix_desc_t *spmtx, const SolverCblk *cblk, int side)
TODO.
static parsec_data_t * parsec_sparse_matrix_data_of(parsec_data_collection_t *mat,...)
Return the data handler associated to the piece of data (uplo, cblknum, bloknum).
static int32_t parsec_sparse_matrix_vpid_of(parsec_data_collection_t *mat,...)
Return the rank of the virtual process owner of the piece of data (uplo, cblknum, bloknum).
static uint32_t parsec_sparse_matrix_rank_of_key(parsec_data_collection_t *mat, parsec_data_key_t key)
Return the rank of the owner of the piece of data (key)
void parsec_sparse_matrix_destroy(parsec_sparse_matrix_desc_t *desc)
Free the PaRSEC descriptor of the sparse matrix.
static void spm_data_key_to_value(parsec_data_key_t key, const SolverMatrix *solvmtx, int *uplo, pastix_int_t *cblknum, pastix_int_t *bloknum)
Compute the triplet (uplo, cblknum, bloknum) from the key.
PaRSEC descriptor stucture for the sparse matrix.
Definition: pastix_parsec.h:35
pastix_int_t cblkmin2d
Definition: solver.h:215
static pastix_int_t blok_rownbr(const SolverBlok *blok)
Compute the number of rows of a block.
Definition: solver.h:390
void * ucoeftab
Definition: solver.h:172
pastix_int_t lrownum
Definition: solver.h:143
void * handler[2]
Definition: solver.h:138
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 fcblknm
Definition: solver.h:140
pastix_int_t cblknbr
Definition: solver.h:208
pastix_int_t cblkmaxblk
Definition: solver.h:216
SolverBlok *restrict bloktab
Definition: solver.h:223
pastix_int_t frownum
Definition: solver.h:142
pastix_int_t coefind
Definition: solver.h:144
SolverBlok * fblokptr
Definition: solver.h:163
pastix_int_t bloknbr
Definition: solver.h:220
pastix_lrblock_t * LRblock[2]
Definition: solver.h:150
SolverCblk *restrict cblktab
Definition: solver.h:222
pastix_int_t stride
Definition: solver.h:164
void * handler[2]
Definition: solver.h:173
int8_t cblktype
Definition: solver.h:159
int ownerid
Definition: solver.h:175
void * lcoeftab
Definition: solver.h:171
Solver block structure.
Definition: solver.h:137
Solver column block structure.
Definition: solver.h:156
Solver column block structure.
Definition: solver.h:200