Go to the documentation of this file.
34 #if defined(PASTIX_WITH_PARSEC)
38 #if defined(PASTIX_WITH_STARPU)
50 struct coeftabinit_s {
51 const SolverMatrix *datacode;
79 struct coeftabinit_s *ciargs = (
struct coeftabinit_s*)args;
80 const SolverMatrix *datacode = ciargs->datacode;
82 const char *dirname = ciargs->dirname;
84 pastix_int_t i, itercblk;
107 for (i=0; i < datacode->ttsknbr[rank]; i++)
109 task = datacode->ttsktab[rank][i];
110 itercblk = datacode->tasktab[task].cblknum;
113 initfunc( side, datacode, bcsc, itercblk, dirname );
141 struct coeftabinit_s args;
143 args.datacode = pastix_data->solvmatr;
144 args.bcsc = pastix_data->bcsc;
147 #if defined(PASTIX_DEBUG_DUMP_COEFTAB)
151 args.dirname = pastix_data->dir_local;
155 SolverCblk *cblk = pastix_data->solvmatr->cblktab;
156 SolverBlok *blok = pastix_data->solvmatr->bloktab;
159 for (i=0; i<pastix_data->solvmatr->cblknbr; i++, cblk++) {
163 for (i=0; i<pastix_data->solvmatr->bloknbr; i++, blok++) {
167 isched_parallel_call( pastix_data->isched,
pcoeftabInit, &args );
189 #if defined(PASTIX_WITH_PARSEC)
191 if ( solvmtx->parsec_desc != NULL ) {
193 free( solvmtx->parsec_desc );
195 solvmtx->parsec_desc = NULL;
198 #if defined(PASTIX_WITH_STARPU)
200 if ( solvmtx->starpu_desc != NULL ) {
202 free( solvmtx->starpu_desc );
204 solvmtx->starpu_desc = NULL;
209 if( solvmtx->cblktab )
213 for ( i = 0; i < solvmtx->cblknbr; i++, cblk++ ) {
215 if( cblk->
cblktype & (CBLK_FANIN|CBLK_RECV) ){
226 struct coeftabcomp_s {
227 SolverMatrix *solvmtx;
229 pastix_atomic_lock_t lock;
255 struct coeftabcomp_s *ccargs = (
struct coeftabcomp_s*)args;
256 SolverMatrix *solvmtx = ccargs->solvmtx;
258 pastix_atomic_lock_t *lock = &(ccargs->lock);
259 pastix_int_t *fullgain = &(ccargs->gain);
262 pastix_int_t i, itercblk;
263 pastix_int_t task, gain = 0;
264 int rank = ctx->rank;
269 case PastixComplex32:
272 case PastixComplex64:
284 for (i=0; i < solvmtx->ttsknbr[rank]; i++)
286 task = solvmtx->ttsktab[rank][i];
287 itercblk = solvmtx->tasktab[task].cblknum;
288 cblk = solvmtx->cblktab + itercblk;
290 if ( cblk->
cblktype & CBLK_COMPRESSED ) {
291 gain += compfunc( solvmtx, side, -1, cblk );
295 pastix_atomic_lock( lock );
297 pastix_atomic_unlock( lock );
322 struct coeftabcomp_s args;
325 args.solvmtx = pastix_data->solvmatr;
326 args.flttype = pastix_data->bcsc->flttype;
327 args.lock = PASTIX_ATOMIC_UNLOCKED;
331 lr = &(pastix_data->solvmatr->lowrank);
337 isched_parallel_call( pastix_data->isched,
pcoeftabComp, (
void*)(&args) );
342 #if defined(PASTIX_WITH_MPI)
364 coeftab_scatter( SolverMatrix *solvmtx,
376 size_t eltsize = pastix_size_of( flttype );
384 cblk = solvmtx->cblktab;
385 for(i=0; i<solvmtx->cblknbr; i++, cblk++)
390 if ( (solvmtx->clustnum == root) &&
391 (solvmtx->clustnum != cblk->
ownerid) )
393 MPI_Send( cblk->
lcoeftab, cblksize, MPI_CHAR,
395 free_fct( side, cblk );
399 if ( (solvmtx->clustnum != root) &&
400 (solvmtx->clustnum == cblk->
ownerid) )
402 MPI_Recv( cblk->
lcoeftab, cblksize, MPI_CHAR,
403 root, i, comm, &status );
429 coeftab_gather( SolverMatrix *solvmtx,
441 size_t eltsize = pastix_size_of( flttype );
444 void *recvbuf = NULL;
451 cblk = solvmtx->cblktab;
452 for( i=0; i<solvmtx->cblknbr; i++, cblk++ )
456 if ( (solvmtx->clustnum == root) &&
457 (solvmtx->clustnum != cblk->
ownerid) )
459 if ( cblk->
cblktype & CBLK_COMPRESSED ) {
465 MALLOC_INTERN( recvbuf, cblksize,
char );
466 MPI_Recv( recvbuf, cblksize, MPI_CHAR,
467 cblk->
ownerid, i, comm, &status );
470 case PastixComplex64:
473 case PastixComplex32:
488 alloc_fct( side, cblk );
489 MPI_Recv( cblk->
lcoeftab, cblksize, MPI_CHAR,
490 cblk->
ownerid, i, comm, &status );
494 if ( (solvmtx->clustnum != root) &&
495 (solvmtx->clustnum == cblk->
ownerid) )
497 if ( cblk->
cblktype & CBLK_COMPRESSED ) {
501 case PastixComplex64:
506 case PastixComplex32:
525 MPI_Send( buffer, bufsize, MPI_CHAR,
529 MPI_Send( cblk->
lcoeftab, cblksize, MPI_CHAR,
548 coeftab_nullify( SolverMatrix *solvmtx )
554 cblk = solvmtx->cblktab;
555 for( i=0; i < solvmtx->cblknbr; i++, cblk++ )
557 if ( solvmtx->clustnum != cblk->
ownerid )
void coeftab_cmemory(SolverMatrix *solvmtx, pastix_fixdbl_t *dparm)
Compute the memory gain of the low-rank form over the full-rank form for the entire matrix.
static pastix_int_t cblk_colnbr(const SolverCblk *cblk)
Compute the number of columns in a column block.
size_t cpucblk_ccompute_size(pastix_coefside_t side, const SolverCblk *cblk)
Compute the size of the buffer to send.
Structure to define the type of function to use for the low-rank kernels and their parameters.
pastix_int_t cpucblk_dcompress(const SolverMatrix *solvmtx, pastix_coefside_t side, int max_ilulvl, SolverCblk *cblk)
Compress a single column block from full-rank to low-rank format.
void cpucblk_sfree(pastix_coefside_t side, SolverCblk *cblk)
Free the cblk structure that store the coefficient.
@ DPARM_COMPRESS_TOLERANCE
void cpucblk_zinit(pastix_coefside_t side, const SolverMatrix *solvmtx, const pastix_bcsc_t *bcsc, pastix_int_t itercblk, const char *directory)
Fully initialize a single cblk.
void(* coeftab_fct_memory_t)(SolverMatrix *, pastix_fixdbl_t *)
Type of the memory gain functions.
Internal column block distributed CSC matrix.
void coeftab_zmemory(SolverMatrix *solvmtx, pastix_fixdbl_t *dparm)
Compute the memory gain of the low-rank form over the full-rank form for the entire matrix.
void cpucblk_sunpack_lr(pastix_coefside_t side, SolverCblk *cblk, void *buffer)
Unpack low rank data and fill the column block concerned by the computation.
Solver column block structure.
enum pastix_coefside_e pastix_coefside_t
Data blocks used in the kernel.
void cpucblk_dfree(pastix_coefside_t side, SolverCblk *cblk)
Free the cblk structure that store the coefficient.
size_t cpucblk_scompute_size(pastix_coefside_t side, const SolverCblk *cblk)
Compute the size of the buffer to send.
void * cpucblk_spack_lr(pastix_coefside_t side, SolverCblk *cblk, size_t size)
Pack low-rank data for column block.
void pastix_gendirectories(pastix_data_t *pastix_data)
Generate a unique temporary directory to store output files.
void cpucblk_calloc(pastix_coefside_t side, SolverCblk *cblk)
Allocate the cblk structure to store the coefficient.
void cpucblk_dinit(pastix_coefside_t side, const SolverMatrix *solvmtx, const pastix_bcsc_t *bcsc, pastix_int_t itercblk, const char *directory)
Fully initialize a single cblk.
pastix_int_t coeftabCompress(pastix_data_t *pastix_data)
Compress the factorized matrix structure if not already done.
static void pcoeftabComp(isched_thread_t *ctx, void *args)
Internal routine called by each static thread to Initialize the solver matrix structure.
@ IPARM_COMPRESS_MIN_HEIGHT
void * cpucblk_dpack_lr(pastix_coefside_t side, SolverCblk *cblk, size_t size)
Pack low-rank data for column block.
size_t cpucblk_zcompute_size(pastix_coefside_t side, const SolverCblk *cblk)
Compute the size of the buffer to send.
void cpucblk_dalloc(pastix_coefside_t side, SolverCblk *cblk)
Allocate the cblk structure to store the coefficient.
void parsec_sparse_matrix_destroy(parsec_sparse_matrix_desc_t *desc)
Free the PaRSEC descriptor of the sparse matrix.
@ IPARM_COMPRESS_MIN_WIDTH
size_t cpucblk_dcompute_size(pastix_coefside_t side, const SolverCblk *cblk)
Compute the size of the buffer to send.
pastix_int_t cpucblk_scompress(const SolverMatrix *solvmtx, pastix_coefside_t side, int max_ilulvl, SolverCblk *cblk)
Compress a single column block from full-rank to low-rank format.
void coeftab_dmemory(SolverMatrix *solvmtx, pastix_fixdbl_t *dparm)
Compute the memory gain of the low-rank form over the full-rank form for the entire matrix.
pastix_int_t cpucblk_ccompress(const SolverMatrix *solvmtx, pastix_coefside_t side, int max_ilulvl, SolverCblk *cblk)
Compress a single column block from full-rank to low-rank format.
void cpucblk_dunpack_lr(pastix_coefside_t side, SolverCblk *cblk, void *buffer)
Unpack low rank data and fill the column block concerned by the computation.
void cpucblk_salloc(pastix_coefside_t side, SolverCblk *cblk)
Allocate the cblk structure to store the coefficient.
void cpucblk_sinit(pastix_coefside_t side, const SolverMatrix *solvmtx, const pastix_bcsc_t *bcsc, pastix_int_t itercblk, const char *directory)
Fully initialize a single cblk.
void cpucblk_zalloc(pastix_coefside_t side, SolverCblk *cblk)
Allocate the cblk structure to store the coefficient.
void cpucblk_cinit(pastix_coefside_t side, const SolverMatrix *solvmtx, const pastix_bcsc_t *bcsc, pastix_int_t itercblk, const char *directory)
Fully initialize a single cblk.
void coeftabExit(SolverMatrix *solvmtx)
Free the solver matrix structure.
void coeftab_smemory(SolverMatrix *solvmtx, pastix_fixdbl_t *dparm)
Compute the memory gain of the low-rank form over the full-rank form for the entire matrix.
void * cpucblk_zpack_lr(pastix_coefside_t side, SolverCblk *cblk, size_t size)
Pack low-rank data for column block.
void pcoeftabInit(isched_thread_t *ctx, void *args)
Internal routine called by each static thread to Initialize the solver matrix structure.
void * cpucblk_cpack_lr(pastix_coefside_t side, SolverCblk *cblk, size_t size)
Pack low-rank data for column block.
coeftab_fct_memory_t coeftabMemory[4]
List of functions to compute the memory gain in low-rank per precision.
void cpucblk_cunpack_lr(pastix_coefside_t side, SolverCblk *cblk, void *buffer)
Unpack low rank data and fill the column block concerned by the computation.
pastix_int_t cpucblk_zcompress(const SolverMatrix *solvmtx, pastix_coefside_t side, int max_ilulvl, SolverCblk *cblk)
Compress a single column block from full-rank to low-rank format.
void coeftabInit(pastix_data_t *pastix_data, pastix_coefside_t side)
Initialize the solver matrix structure.
#define pastix_coeftype_t
Arithmetic types.
void cpucblk_cfree(pastix_coefside_t side, SolverCblk *cblk)
Free the cblk structure that store the coefficient.
void cpucblk_zfree(pastix_coefside_t side, SolverCblk *cblk)
Free the cblk structure that store the coefficient.
void cpucblk_zunpack_lr(pastix_coefside_t side, SolverCblk *cblk, void *buffer)
Unpack low rank data and fill the column block concerned by the computation.
void starpu_sparse_matrix_destroy(starpu_sparse_matrix_desc_t *desc)
Free the StarPU descriptor of the sparse matrix.