20 #include "common/common.h"
29 #if defined( PASTIX_WITH_MPI )
52 SolverMatrix *solvmtx,
60 assert( cblk->
cblktype & CBLK_FANIN );
65 #if defined(PASTIX_DEBUG_MPI)
66 fprintf( stderr,
"[%2d] Post Isend for cblk %ld toward %2d ( %ld Bytes )\n",
70 rc = MPI_Isend( buffer, bufsize, MPI_CHAR,
72 assert( rc == MPI_SUCCESS );
74 solverCommMatrixAdd( solvmtx, cblk->
ownerid, bufsize );
77 pastix_atomic_lock( &(solvmtx->reqlock) );
79 assert( solvmtx->reqidx[ solvmtx->reqnum ] == -1 );
80 assert( solvmtx->reqnum >= 0 );
81 assert( solvmtx->reqnum < solvmtx->reqnbr );
83 solvmtx->reqtab[ solvmtx->reqnum ] = request;
84 solvmtx->reqidx[ solvmtx->reqnum ] = cblk - solvmtx->cblktab;
87 pastix_atomic_unlock( &(solvmtx->reqlock) );
114 const SolverMatrix *solvmtx,
117 assert( cblk->
cblktype & CBLK_FANIN );
119 #if defined(PASTIX_DEBUG_MPI)
122 cblksize = cblksize *
cblk_colnbr( cblk ) * cblk->
stride *
sizeof(pastix_complex32_t);
124 fprintf( stderr,
"[%2d] Isend for cblk %ld toward %2d ( %ld Bytes ) (DONE)\n",
125 solvmtx->clustnum, (
long)cblk->
gcblknum, cblk->
ownerid, (
long)cblksize );
155 SolverMatrix *solvmtx,
157 const MPI_Status *status,
161 int src = status->MPI_SOURCE;
162 int tag = status->MPI_TAG;
164 assert( ( 0 <= src ) && ( src < solvmtx->clustnbr ) );
165 assert( ( 0 <= tag ) && ( tag < solvmtx->gcblknbr ) );
170 fcbk = solvmtx->cblktab + solvmtx->gcbl2loc[ tag ];
174 while( cblk->
ownerid != src ) {
176 assert( cblk >= solvmtx->cblktab );
178 assert( cblk->
cblktype & CBLK_RECV );
181 #if defined(PASTIX_DEBUG_MPI)
183 pastix_int_t size = (
cblk_colnbr(cblk) * cblk->
stride) *
sizeof(pastix_complex32_t);
190 MPI_Get_count( status, MPI_CHAR, &count );
191 assert( (cblk->
cblktype & CBLK_COMPRESSED) ||
192 (!(cblk->
cblktype & CBLK_COMPRESSED) && (count == size)) );
195 fprintf( stderr,
"[%2d] Irecv of size %d/%ld for cblk %ld (DONE)\n",
196 solvmtx->clustnum, count, (
long)size, (
long)cblk->
gcblknum );
253 SolverMatrix *solvmtx,
257 const MPI_Status *statuses )
259 pastix_int_t i, reqid;
260 int nbrequest = outcount;
262 for( i = 0; i < outcount; i++ ){
268 if ( solvmtx->reqidx[reqid] == -1 ) {
274 memcpy( &status, statuses + i,
sizeof(MPI_Status) );
275 MPI_Get_count( &status, MPI_CHAR, &size );
277 MALLOC_INTERN( recvbuf, size,
char );
278 memcpy( recvbuf, solvmtx->rcoeftab, size );
283 if ( solvmtx->recvcnt > 0 ) {
284 MPI_Start( solvmtx->reqtab + reqid );
288 MPI_Request_free( solvmtx->reqtab + reqid );
289 solvmtx->reqtab[reqid] = MPI_REQUEST_NULL;
292 cpucblk_crequest_handle_recv( side, solvmtx, threadid,
299 SolverCblk *cblk = solvmtx->cblktab + solvmtx->reqidx[ reqid ];
300 assert( cblk->
cblktype & CBLK_FANIN );
302 cpucblk_crequest_handle_fanin( side, solvmtx, cblk );
305 solvmtx->reqidx[ reqid ] = -1;
326 cpucblk_cupdate_reqtab( SolverMatrix *solvmtx )
329 MPI_Request *outrequest = solvmtx->reqtab;
330 pastix_int_t *outreqloc = solvmtx->reqidx;
334 MPI_Request *inrequest = solvmtx->reqtab;
335 pastix_int_t *inreqloc = solvmtx->reqidx;
339 while( (outreqnbr < solvmtx->reqnum) &&
340 (*outrequest != MPI_REQUEST_NULL) )
347 inrequest = outrequest;
348 inreqloc = outreqloc;
349 inreqnbr = outreqnbr;
350 for( ; inreqnbr < solvmtx->reqnum;
351 inreqnbr++, inrequest++, inreqloc++ )
353 if ( *inrequest == MPI_REQUEST_NULL )
359 *outrequest = *inrequest;
360 *outreqloc = *inreqloc;
370 memset( outreqloc, 0xff, (solvmtx->reqnbr - outreqnbr) *
sizeof(pastix_int_t) );
373 #if defined(PASTIX_DEBUG_MPI)
375 for( i = outreqnbr; i < solvmtx->reqnbr; i++ )
377 solvmtx->reqtab[i] = MPI_REQUEST_NULL;
380 assert( outreqnbr < solvmtx->reqnum );
381 solvmtx->reqnum = outreqnbr;
410 SolverMatrix *solvmtx,
413 pthread_t tid = pthread_self();
415 int nbrequest, nbfree;
416 int indexes[ solvmtx->reqnbr ];
417 MPI_Status statuses[ solvmtx->reqnbr ];
420 pthread_mutex_lock( &pastix_comm_lock );
421 if ( pastix_comm_tid == (pthread_t)-1 ) {
422 pastix_comm_tid = tid;
424 pthread_mutex_unlock( &pastix_comm_lock );
426 if ( tid != pastix_comm_tid ) {
439 pastix_atomic_lock( &(solvmtx->reqlock) );
440 nbrequest = solvmtx->reqnum;
441 pastix_atomic_unlock( &(solvmtx->reqlock) );
443 while( (outcount > 0) && (nbrequest > 0) )
445 MPI_Testsome( nbrequest, solvmtx->reqtab, &outcount, indexes, statuses );
449 if ( outcount > 0 ) {
450 nbfree = cpucblk_crequest_handle( side, solvmtx, threadid,
451 outcount, indexes, statuses );
458 pastix_atomic_lock( &(solvmtx->reqlock) );
460 cpucblk_cupdate_reqtab( solvmtx );
462 nbrequest = solvmtx->reqnum;
463 pastix_atomic_unlock( &(solvmtx->reqlock) );
466 pastix_comm_tid = -1;
501 SolverMatrix *solvmtx,
504 #if defined(PASTIX_WITH_MPI)
505 if ( cblk->
cblktype & CBLK_FANIN ) {
520 cpucblk_cmpi_progress( side, solvmtx, rank );
523 assert( !(cblk->
cblktype & (CBLK_FANIN | CBLK_RECV)) );
524 do { }
while( cblk->
ctrbcnt > 0 );
559 SolverMatrix *solvmtx,
564 ctrbcnt = pastix_atomic_dec_32b( &(fcbk->
ctrbcnt) );
566 #if defined(PASTIX_WITH_MPI)
567 if ( fcbk->
cblktype & CBLK_FANIN ) {
568 cpucblk_cisend( side, solvmtx, fcbk );
574 if ( solvmtx->computeQueue ) {
611 SolverMatrix *solvmtx )
619 #if defined(PASTIX_WITH_MPI)
623 int reqnbr = solvmtx->reqnum;
626 #if defined(PASTIX_DEBUG_MPI)
627 fprintf( stderr,
"[%2d] Wait for all pending communications\n",
631 for( i=0; i<reqnbr; i++ )
633 if ( solvmtx->reqtab[i] == MPI_REQUEST_NULL ) {
640 assert( solvmtx->reqidx[i] != -1 );
642 rc = MPI_Wait( solvmtx->reqtab + i, &status );
643 assert( rc == MPI_SUCCESS );
645 cblk = solvmtx->cblktab + solvmtx->reqidx[i];
648 assert( cblk->
cblktype & CBLK_FANIN );
650 cpucblk_crequest_handle_fanin( side, solvmtx, cblk );
654 assert( solvmtx->reqnum == 0 );