PaStiX Handbook  6.2.1
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-2021 Bordeaux INP, CNRS (LaBRI UMR 5800), Inria,
8  * Univ. Bordeaux. All rights reserved.
9  *
10  * @version 6.2.0
11  * @author Mathieu Faverge
12  * @author Pierre Ramet
13  * @date 2021-01-03
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 void
448 parsec_sparse_matrix_init_fermi( parsec_sparse_matrix_desc_t *spmtx,
449  const SolverMatrix *solvmtx )
450 {
451  gpu_device_t* gpu_device;
452  SolverBlok *blok;
453  pastix_int_t i, b, bloknbr, ndevices;
454  size_t size;
455  int *tmp, *bloktab;
456 
457  ndevices = parsec_devices_enabled();
458  if ( ndevices <= 2 )
459  return;
460 
461  bloknbr = solvmtx->bloknbr;
462  size = 2 * bloknbr * sizeof(int);
463 
464  /**
465  * Initialize array on CPU
466  */
467  bloktab = (int*)malloc( size );
468  tmp = bloktab;
469  for (b=0, blok = solvmtx->bloktab;
470  b < bloknbr;
471  b++, blok++, tmp+=2)
472  {
473  tmp[0] = blok->frownum;
474  tmp[1] = blok->lrownum;
475  }
476 
477  ndevices -= 2;
478  spmtx->gpu_blocktab = calloc(ndevices, sizeof(void*));
479 
480  fprintf(stderr, "ndevices = %ld\n", ndevices );
481  for(i = 0; i < ndevices; i++) {
482  if( NULL == (gpu_device = (gpu_device_t*)parsec_devices_get(i+2)) ) continue;
483 
484  fprintf(stderr, "cuda index = %d\n", gpu_device->cuda_index );
485  cudaSetDevice( gpu_device->cuda_index );
486 
487  cudaMalloc( &(spmtx->gpu_blocktab[i]),
488  size );
489 
490  cudaMemcpy( spmtx->gpu_blocktab[i],
491  bloktab, size,
492  cudaMemcpyHostToDevice );
493  }
494  free(bloktab);
495 }
496 
497 void
498 parsec_sparse_matrix_destroy_fermi( parsec_sparse_matrix_desc_t *spmtx )
499 {
500  gpu_device_t* gpu_device;
501  pastix_int_t i, ndevices;
502 
503  ndevices = parsec_devices_enabled();
504  if ( ndevices <= 2 )
505  return;
506 
507  ndevices -= 2;
508  for(i = 0; i < ndevices; i++) {
509  if( NULL == (gpu_device = (gpu_device_t*)parsec_devices_get(i+2)) ) continue;
510 
511  cudaSetDevice( gpu_device->cuda_index );
512  cudaFree( spmtx->gpu_blocktab[i] );
513  }
514 
515  free( spmtx->gpu_blocktab );
516 }
517 #endif /*defined(PASTIX_CUDA_FERMI)*/
518 
519 static inline void
520 pastix_parsec_register_cblk_lr( parsec_data_collection_t *o,
521  parsec_data_t **handler,
522  pastix_int_t id,
523  const SolverCblk *cblk,
524  int side )
525 {
526  pastix_lrblock_t *dataptr = cblk->fblokptr->LRblock[side];
527  size_t size = ( cblk[1].fblokptr - cblk[0].fblokptr ) * sizeof( pastix_lrblock_t );
528 
529  parsec_data_create( handler, o, id, dataptr, size );
530 }
531 
532 static inline void
533 pastix_parsec_register_cblk_fr( parsec_data_collection_t *o,
534  parsec_data_t **handler,
535  pastix_int_t id,
536  const parsec_sparse_matrix_desc_t *spmtx,
537  const SolverCblk *cblk,
538  int side )
539 {
540  void *dataptr = ( side == PastixUCoef ) ? cblk->ucoeftab : cblk->lcoeftab;
541  size_t size = (size_t)cblk->stride * (size_t)cblk_colnbr( cblk ) * (size_t)spmtx->typesze;
542 
543  parsec_data_create( handler, o, id, dataptr, size );
544 }
545 
546 static inline void
547 pastix_parsec_register_cblk( parsec_data_collection_t *o,
548  pastix_int_t cblknum,
549  const parsec_sparse_matrix_desc_t *spmtx,
550  const SolverCblk *cblk )
551 {
552  parsec_data_t **handler = (parsec_data_t **)( cblk->handler );
553 
554  if ( cblk->cblktype & CBLK_COMPRESSED ) {
555  pastix_parsec_register_cblk_lr( o, handler, cblknum * 2, cblk, PastixLCoef );
556 
557  if ( spmtx->mtxtype == PastixGeneral ) {
558  pastix_parsec_register_cblk_lr( o, handler + 1, cblknum * 2 + 1, cblk, PastixUCoef );
559  }
560  }
561  else {
562  pastix_parsec_register_cblk_fr( o, handler, cblknum * 2, spmtx, cblk, PastixLCoef );
563  if ( spmtx->mtxtype == PastixGeneral ) {
564  pastix_parsec_register_cblk_fr(
565  o, handler + 1, cblknum * 2 + 1, spmtx, cblk, PastixUCoef );
566  }
567  }
568 
569  assert(cblk->handler[0]);
570  assert(cblk->handler[1] || ( spmtx->mtxtype != PastixGeneral ));
571 }
572 
573 /**
574  *******************************************************************************
575  *
576  * @brief Generate the PaRSEC descriptor of the sparse matrix.
577  *
578  * This function creates the PaRSEC descriptor that will provide tha data
579  * mapping and memory location to PaRSEC for the computation.
580  *
581  *******************************************************************************
582  *
583  * @param[inout] solvmtx
584  * The solver matrix structure that describes the sparse matrix for
585  * PaStiX.
586  *
587  * @param[in] typesize
588  * The memory size of the arithmetic used to store the matrix
589  * coefficients.
590  *
591  * @param[in] mtxtype
592  * The type of sparse matrix to describe.
593  * @arg PastixGeneral: The sparse matrix is general.
594  * @arg PastixSymmetric: The sparse matrix is lower triangular symmetric.
595  * @arg PastixHermitian: The sparse matrix is lower triangular hermitian.
596  *
597  * @param[in] nodes
598  * The number of processes used to solve the problem.
599  *
600  * @param[in] myrank
601  * The rank of the calling process.
602  *
603  ******************************************************************************/
604 void
605 parsec_sparse_matrix_init( SolverMatrix *solvmtx,
606  int typesize, int mtxtype,
607  int nodes, int myrank )
608 {
609  parsec_sparse_matrix_desc_t *spmtx = solvmtx->parsec_desc;
610  parsec_data_collection_t *o;
611  pastix_int_t cblknbr, cblkmin2d, ld;
612  parsec_data_key_t key1, key2;
613  SolverCblk *cblk;
614  SolverBlok *blok, *fblok, *lblok;
615  pastix_int_t m=0, n=0, cblknum, nbrow;
616  size_t size, offset;
617  char *ptrL, *ptrU;
618  pastix_lrblock_t *dataptrL, *dataptrU;
619 
620  if ( spmtx != NULL ) {
622  }
623  else {
625  }
626 
627  o = (parsec_data_collection_t*)spmtx;
628  parsec_data_collection_init( o, nodes, myrank );
629 
630  o->data_key = parsec_sparse_matrix_data_key;
631 #if defined(PARSEC_PROF_TRACE)
632  o->key_to_string = parsec_sparse_matrix_key_to_string;
633 #endif
634 
635  o->rank_of = parsec_sparse_matrix_rank_of;
636  o->rank_of_key = parsec_sparse_matrix_rank_of_key;
637  o->vpid_of = parsec_sparse_matrix_vpid_of;
638  o->vpid_of_key = parsec_sparse_matrix_vpid_of_key;
639  o->data_of = parsec_sparse_matrix_data_of;
640  o->data_of_key = parsec_sparse_matrix_data_of_key;
641 
642  spmtx->typesze = typesize;
643  spmtx->mtxtype = mtxtype;
644  spmtx->solvmtx = solvmtx;
645 
646  cblknbr = solvmtx->cblknbr;
647  cblkmin2d = solvmtx->cblkmin2d;
648  ld = solvmtx->cblkmaxblk * 2;
649  key1 = 2 * cblknbr;
650 
651  /* Initialize 1D cblk handlers */
652  cblk = spmtx->solvmtx->cblktab;
653  for(cblknum = 0;
654  cblknum < cblkmin2d;
655  cblknum++, n++, cblk++ )
656  {
657  if ( cblk->ownerid != myrank ) {
658  continue;
659  } else {
660  pastix_parsec_register_cblk( o, cblknum, spmtx, cblk );
661  }
662  }
663 
664  /* Initialize 2D cblk handlers */
665  cblk = spmtx->solvmtx->cblktab + cblkmin2d;
666  for(cblknum = cblkmin2d, n = 0;
667  cblknum < cblknbr;
668  cblknum++, n++, cblk++ )
669  {
670  if ( cblk->ownerid != myrank ) {
671  continue;
672  }
673 
674  pastix_parsec_register_cblk( o, cblknum, spmtx, cblk );
675 
676  if ( !(cblk->cblktype & CBLK_TASKS_2D) ) {
677  continue;
678  }
679 
680  if ( cblk->cblktype & CBLK_COMPRESSED ) {
681  /*
682  * Diagonal block
683  */
684  dataptrL = cblk->fblokptr->LRblock[0];
685  dataptrU = cblk->fblokptr->LRblock[1];
686  blok = cblk->fblokptr;
687  size = sizeof( pastix_lrblock_t );
688  offset = 0;
689  key2 = n * ld;
690 
691  assert( offset == 0 );
692  parsec_data_create( (parsec_data_t **)&( blok->handler[0] ),
693  o, key1 + key2, dataptrL + offset, size );
694 
695  if ( mtxtype == PastixGeneral ) {
696  parsec_data_create( (parsec_data_t **)&( blok->handler[1] ),
697  o, key1 + key2 + 1, dataptrU + offset, size );
698  }
699  else {
700  blok->handler[1] = NULL;
701  }
702 
703  /*
704  * Off-diagonal blocks
705  */
706  blok++;
707  key2 += 2;
708  lblok = cblk[1].fblokptr;
709  for ( ; blok < lblok; blok++, key2 += 2 ) {
710  fblok = blok;
711  m = 0;
712  size = blok_rownbr( blok );
713  offset = blok - cblk->fblokptr;
714  nbrow = 1;
715 
716  while ( ( blok < lblok ) && ( blok[0].fcblknm == blok[1].fcblknm ) &&
717  ( blok[0].lcblknm == blok[1].lcblknm ) ) {
718  blok++;
719  m++;
720  nbrow++;
721  }
722  size = nbrow;
723 
724  parsec_data_create( (parsec_data_t **)&( fblok->handler[0] ),
725  &spmtx->super,
726  key1 + key2,
727  dataptrL + offset,
728  size );
729 
730  if ( mtxtype == PastixGeneral ) {
731  parsec_data_create( (parsec_data_t **)&( fblok->handler[1] ),
732  &spmtx->super,
733  key1 + key2 + 1,
734  dataptrU + offset,
735  size );
736  }
737  else {
738  fblok->handler[1] = NULL;
739  }
740 
741  key2 += m * 2;
742  }
743  }
744  else {
745  /*
746  * Diagonal block
747  */
748  ptrL = cblk->lcoeftab;
749  ptrU = cblk->ucoeftab;
750  blok = cblk->fblokptr;
751  size = blok_rownbr( blok ) * cblk_colnbr( cblk ) * (size_t)spmtx->typesze;
752  offset = blok->coefind * (size_t)spmtx->typesze;
753  key2 = n * ld;
754 
755  assert( offset == 0 );
756  parsec_data_create(
757  (parsec_data_t **)&( blok->handler[0] ), o, key1 + key2, ptrL + offset, size );
758 
759  if ( mtxtype == PastixGeneral ) {
760  parsec_data_create( (parsec_data_t **)&( blok->handler[1] ),
761  o,
762  key1 + key2 + 1,
763  ptrU + offset,
764  size );
765  }
766  else {
767  blok->handler[1] = NULL;
768  }
769 
770  /*
771  * Off-diagonal blocks
772  */
773  blok++;
774  key2 += 2;
775  lblok = cblk[1].fblokptr;
776  for ( ; blok < lblok; blok++, key2 += 2 ) {
777  fblok = blok;
778  m = 0;
779  size = blok_rownbr( blok );
780  offset = blok->coefind * (size_t)spmtx->typesze;
781 
782  while ( ( blok < lblok ) &&
783  ( blok[0].fcblknm == blok[1].fcblknm ) &&
784  ( blok[0].lcblknm == blok[1].lcblknm ) )
785  {
786  blok++;
787  m++;
788  size += blok_rownbr( blok );
789  }
790  size *= cblk_colnbr( cblk ) * (size_t)spmtx->typesze;
791 
792  parsec_data_create( (parsec_data_t **)&( fblok->handler[0] ),
793  &spmtx->super,
794  key1 + key2,
795  ptrL + offset,
796  size );
797 
798  if ( mtxtype == PastixGeneral ) {
799  parsec_data_create( (parsec_data_t **)&( fblok->handler[1] ),
800  &spmtx->super,
801  key1 + key2 + 1,
802  ptrU + offset,
803  size );
804  }
805  else {
806  fblok->handler[1] = NULL;
807  }
808 
809  key2 += m * 2;
810  }
811  }
812  }
813 
814 #if defined( PASTIX_CUDA_FERMI )
815  parsec_sparse_matrix_init_fermi( spmtx, solvmtx );
816 #endif
817 
818  solvmtx->parsec_desc = spmtx;
819 }
820 
821 /**
822  *******************************************************************************
823  *
824  * @brief Free the PaRSEC descriptor of the sparse matrix.
825  *
826  * This function destroys the PaRSEC descriptor, but do not free the matrix data
827  * that are managed by PaStiX.
828  *
829  *******************************************************************************
830  *
831  * @param[inout] spmtx
832  * The descriptor to free.
833  *
834  ******************************************************************************/
835 void
837 {
838  SolverCblk *cblk;
839  SolverBlok *blok;
840  pastix_int_t i, cblkmin2d;
841 
842 #if defined(PASTIX_CUDA_FERMI)
843  parsec_sparse_matrix_destroy_fermi( spmtx );
844 #endif
845 
846  cblkmin2d = spmtx->solvmtx->cblkmin2d;
847  cblk = spmtx->solvmtx->cblktab;
848  for(i=0; i<cblkmin2d; i++, cblk++)
849  {
850  if ( cblk->handler[0] ) {
851  parsec_data_destroy( cblk->handler[0] );
852 
853  if ( spmtx->mtxtype == PastixGeneral ) {
854  parsec_data_destroy( cblk->handler[1] );
855  }
856  }
857 
858  cblk->handler[0] = NULL;
859  cblk->handler[1] = NULL;
860  }
861 
862  for(i=cblkmin2d; i<spmtx->solvmtx->cblknbr; i++, cblk++)
863  {
864  if ( cblk->handler[0] ) {
865  parsec_data_destroy( cblk->handler[0] );
866  if ( spmtx->mtxtype == PastixGeneral ) {
867  parsec_data_destroy( cblk->handler[1] );
868  }
869  }
870 
871  cblk->handler[0] = NULL;
872  cblk->handler[1] = NULL;
873 
874  blok = cblk->fblokptr;
875  while( blok < cblk[1].fblokptr )
876  {
877  if ( blok->handler[0] ) {
878  parsec_data_destroy( blok->handler[0] );
879  if ( spmtx->mtxtype == PastixGeneral ) {
880  parsec_data_destroy( blok->handler[1] );
881  }
882  }
883 
884  blok->handler[0] = NULL;
885  blok->handler[1] = NULL;
886 
887  blok++;
888  }
889  }
890 
891  parsec_data_collection_destroy( (parsec_data_collection_t*)spmtx );
892 }
893 
894 /**
895  *@}
896  */
solver_cblk_s::ownerid
int ownerid
Definition: solver.h:146
solver_blok_s::frownum
pastix_int_t frownum
Definition: solver.h:112
solver.h
blok_rownbr
static pastix_int_t blok_rownbr(const SolverBlok *blok)
Compute the number of rows of a block.
Definition: solver.h:313
cblk_colnbr
static pastix_int_t cblk_colnbr(const SolverCblk *cblk)
Compute the number of columns in a column block.
Definition: solver.h:247
solver_cblk_s::fblokptr
SolverBlok * fblokptr
Definition: solver.h:134
parsec_sparse_matrix_desc_s::mtxtype
int mtxtype
Definition: pastix_parsec.h:36
solver_cblk_s::stride
pastix_int_t stride
Definition: solver.h:135
solver_cblk_s
Solver column block structure.
Definition: solver.h:127
parsec_sparse_matrix_rank_of
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).
Definition: parsec_sparse_matrix.c:178
parsec_sparse_matrix_vpid_of_key
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)
Definition: parsec_sparse_matrix.c:285
pastix_parsec.h
solver_blok_s::handler
void * handler[2]
Definition: solver.h:108
solver_blok_s
Solver block structure.
Definition: solver.h:107
pastix_lrblock_s
The block low-rank structure to hold a matrix in low-rank form.
Definition: pastix_lowrank.h:112
PastixGeneral
@ PastixGeneral
Definition: api.h:435
solver_cblk_s::ucoeftab
void * ucoeftab
Definition: solver.h:143
parsec_sparse_matrix_data_of_key
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).
Definition: parsec_sparse_matrix.c:367
pastix_lrblock_t
struct pastix_lrblock_s pastix_lrblock_t
The block low-rank structure to hold a matrix in low-rank form.
parsec_sparse_matrix_destroy
void parsec_sparse_matrix_destroy(parsec_sparse_matrix_desc_t *desc)
Free the PaRSEC descriptor of the sparse matrix.
Definition: parsec_sparse_matrix.c:836
solver_cblk_s::lcoeftab
void * lcoeftab
Definition: solver.h:142
PastixUCoef
@ PastixUCoef
Definition: api.h:457
parsec_sparse_matrix_desc_s::gpu_blocktab
void ** gpu_blocktab
Definition: pastix_parsec.h:38
parsec_sparse_matrix_rank_of_key
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)
Definition: parsec_sparse_matrix.c:219
parsec_sparse_matrix_vpid_of
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).
Definition: parsec_sparse_matrix.c:260
PastixLCoef
@ PastixLCoef
Definition: api.h:456
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
solver_blok_s::coefind
pastix_int_t coefind
Definition: solver.h:114
solver_blok_s::lrownum
pastix_int_t lrownum
Definition: solver.h:113
parsec_sparse_matrix_desc_s
Definition: pastix_parsec.h:33
solver_cblk_s::cblktype
int8_t cblktype
Definition: solver.h:130
parsec_sparse_matrix_data_key
static uint32_t parsec_sparse_matrix_data_key(parsec_data_collection_t *mat,...)
Compute the unique key from the triplet (uplo, cblknum, bloknum).
Definition: parsec_sparse_matrix.c:122
solver_blok_s::LRblock
pastix_lrblock_t * LRblock[2]
Definition: solver.h:121
parsec_sparse_matrix_desc_s::typesze
int typesze
Definition: pastix_parsec.h:35
parsec_sparse_matrix_desc_s::super
parsec_data_collection_t super
Definition: pastix_parsec.h:34
parsec_sparse_matrix_desc_s::solvmtx
SolverMatrix * solvmtx
Definition: pastix_parsec.h:37
spm_data_key_to_value
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.
Definition: parsec_sparse_matrix.c:61
solver_cblk_s::handler
void * handler[2]
Definition: solver.h:144
parsec_sparse_matrix_data_of
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).
Definition: parsec_sparse_matrix.c:316
solver_blok_s::fcblknm
pastix_int_t fcblknm
Definition: solver.h:110