PaStiX Handbook  6.4.0
cpucblk_dmpi_rhs_fwd.c
Go to the documentation of this file.
1 /**
2  *
3  * @file cpucblk_dmpi_rhs_fwd.c
4  *
5  * Precision dependent routines to manag communications for the solve part.
6  *
7  * @copyright 2015-2024 Bordeaux INP, CNRS (LaBRI UMR 5800), Inria,
8  * Univ. Bordeaux. All rights reserved.
9  *
10  * @version 6.4.0
11  * @author Pierre Ramet
12  * @author Mathieu Faverge
13  * @author Tony Delarue
14  * @author Alycia Lisito
15  * @date 2024-07-05
16  *
17  * @generated from /builds/solverstack/pastix/kernels/cpucblk_zmpi_rhs_fwd.c, normal z -> d, Tue Oct 8 14:17:24 2024
18  *
19  **/
20 #include "common/common.h"
21 #include "common/pastixdata.h"
22 #include <lapacke.h>
23 #include "blend/solver.h"
25 #include "kernels.h"
26 #include "pastix_dcores.h"
27 #include "pastix_dlrcores.h"
28 
29 /**
30  *******************************************************************************
31  *
32  * @brief Send the rhs associated to a cblk->lcolidx to the remote node.
33  *
34  *******************************************************************************
35  *
36  * @param[in] solvmtx
37  * The solver matrix holding the communicator.
38  *
39  * @param[in] cblk
40  * The cblk which defines the part to sent.
41  *
42  * @param[in] rhsb
43  * The pointer to the rhs data structure that holds the vectors of the
44  * right hand side.
45  *
46  *******************************************************************************/
47 void
49  SolverCblk *cblk,
50  pastix_rhs_t rhsb )
51 {
52 #if defined(PASTIX_WITH_MPI)
53  double *b;
54  pastix_int_t colnbr = cblk_colnbr(cblk);
55  pastix_int_t size = colnbr * rhsb->n;
56  pastix_int_t idx = - cblk->bcscnum - 1;
57  int rc;
58 
59  assert( colnbr <= solvmtx->colmax );
60  assert( cblk->cblktype & CBLK_FANIN );
61  assert( rhsb->cblkb[ idx ] != NULL );
62 
63 #if defined (PASTIX_DEBUG_MPI)
64  fprintf( stderr, "[%2d] RHS Fwd: Send cblk %ld to %2d at index %ld of size %ld\n",
65  solvmtx->clustnum, (long)cblk->gcblknum, cblk->ownerid,
66  (long)cblk->lcolidx, (long)colnbr );
67 #endif
68 
69  b = (double*)(rhsb->cblkb[ idx ]);
70  assert( b != NULL );
71 
72  rc = MPI_Send( b, size, PASTIX_MPI_DOUBLE,
73  cblk->ownerid, cblk->gcblknum, solvmtx->solv_comm );
74  assert( rc == MPI_SUCCESS );
75 
76  memFree_null( rhsb->cblkb[ idx ] );
77 
78  (void)rc;
79 #else
80  (void)solvmtx;
81  (void)cblk;
82  (void)rhsb;
83 #endif
84 }
85 
86 /**
87  *******************************************************************************
88  *
89  * @brief Send the rhs associated to a cblk->lcolidx to the remote node.
90  *
91  *******************************************************************************
92  *
93  * @param[in] solvmtx
94  * The solver matrix holding the communicator.
95  *
96  * @param[in] cblk
97  * The cblk which defines the part to sent.
98  *
99  * @param[in] rhsb
100  * The pointer to the rhs data structure that holds the vectors of the
101  * right hand side.
102  *
103  *******************************************************************************/
104 void
106  SolverCblk *cblk,
107  pastix_rhs_t rhsb )
108 {
109 #if defined(PASTIX_WITH_MPI)
110  double *b = rhsb->b;
111  pastix_int_t colnbr = cblk_colnbr(cblk);
112  pastix_int_t idx = - cblk->bcscnum - 1;
113  int rc;
114 
115  assert( colnbr <= solvmtx->colmax );
116  assert( cblk->cblktype & CBLK_RECV );
117  assert( rhsb->cblkb[ idx ] == NULL );
118 
119 #if defined (PASTIX_DEBUG_MPI)
120  fprintf( stderr, "[%2d] RHS Bwd: Send cblk %ld to %2d at index %ld of size %ld\n",
121  solvmtx->clustnum, (long)cblk->gcblknum, cblk->ownerid,
122  (long)cblk->lcolidx, (long)colnbr );
123 #endif
124 
125  b += cblk->lcolidx;
126  if ( rhsb->n > 1 ) {
127  rhsb->cblkb[ idx ] = malloc( colnbr * rhsb->n * sizeof(double) );
128  rc = LAPACKE_dlacpy_work( LAPACK_COL_MAJOR, 'A', colnbr, rhsb->n, b, rhsb->ld, rhsb->cblkb[ idx ], colnbr );
129  assert( rc == 0 );
130  b = rhsb->cblkb[ idx ];
131  }
132 
133  rc = MPI_Send( b, colnbr * rhsb->n, PASTIX_MPI_DOUBLE,
134  cblk->ownerid, cblk->gcblknum, solvmtx->solv_comm );
135  assert( rc == MPI_SUCCESS );
136 
137  if ( rhsb->n > 1 ) {
138  memFree_null( rhsb->cblkb[ idx ] );
139  }
140 
141  (void)rc;
142 #else
143  (void)solvmtx;
144  (void)cblk;
145  (void)rhsb;
146 #endif
147 }
148 
149 /**
150  *******************************************************************************
151  *
152  * @brief Receive the rhs associated to a cblk->lcolidx to the remote node.
153  *
154  *******************************************************************************
155  *
156  * @param[in] solvmtx
157  * The solver matrix holding the communicator.
158  *
159  * @param[in] cblk
160  * The cblk which may define the part to sent.
161  *
162  * @param[inout] rhsb
163  * The pointer to the rhs data structure that holds the vectors of the
164  * right hand side.
165  *
166  *******************************************************************************/
167 void
169  SolverCblk *cblk,
170  pastix_rhs_t rhsb )
171 {
172 #if defined(PASTIX_WITH_MPI)
173  MPI_Status status;
174  pastix_int_t colnbr = cblk_colnbr(cblk);
175  pastix_int_t idx = - cblk->bcscnum - 1;
176  int rc;
177 
178  assert( colnbr <= solvmtx->colmax );
179  assert( cblk->cblktype & CBLK_FANIN );
180 
181 #if defined (PASTIX_DEBUG_MPI)
182  fprintf( stderr, "[%2d] RHS Bwd: Recv cblk %ld from %ld at index %ld of size %ld\n",
183  solvmtx->clustnum, (long)cblk->gcblknum, (long)cblk->ownerid,
184  (long)cblk->lcolidx, (long)colnbr );
185 #endif
186 
187  assert( rhsb->cblkb[ idx ] == NULL );
188  rhsb->cblkb[ idx ] = malloc( colnbr * rhsb->n * sizeof(double) );
189 
190  rc = MPI_Recv( rhsb->cblkb[ idx ], colnbr * rhsb->n, PASTIX_MPI_DOUBLE,
191  cblk->ownerid, cblk->gcblknum, solvmtx->solv_comm, &status );
192  assert( rc == MPI_SUCCESS );
193 
194 #if defined (PASTIX_DEBUG_MPI)
195  fprintf( stderr, "[%2d] RHS Bwd: Received cblk %ld from %2d\n",
196  solvmtx->clustnum, (long)cblk->gcblknum, status.MPI_SOURCE );
197 #endif
198 
199  (void)rc;
200 #else
201  (void)solvmtx;
202  (void)cblk;
203  (void)rhsb;
204 #endif
205 }
206 
207 /**
208  *******************************************************************************
209  *
210  * @brief Receive the rhs associated to a cblk->lcolidx to the remote node.
211  *
212  *******************************************************************************
213  *
214  * @param[in] solvmtx
215  * The solver matrix holding the communicator.
216  *
217  * @param[in] cblk
218  * The cblk which may define the part to sent.
219  *
220  * @param[inout] work
221  * The temporary buffer to receive the remote data
222  *
223  * @param[inout] rhsb
224  * The pointer to the rhs data structure that holds the vectors of the
225  * right hand side.
226  *
227  *******************************************************************************/
228 void
230  SolverCblk *cblk,
231  double *work,
232  pastix_rhs_t rhsb )
233 {
234 #if defined(PASTIX_WITH_MPI)
235  double *b = rhsb->b;
236  pastix_int_t colnbr = cblk_colnbr(cblk);
237  MPI_Status status;
238  int rc;
239 
240  assert( colnbr <= solvmtx->colmax );
241  assert( cblk->cblktype & CBLK_RECV );
242 
243 #if defined (PASTIX_DEBUG_MPI)
244  fprintf( stderr, "[%2d] RHS Fwd: Recv cblk %ld from %ld at index %ld of size %ld\n",
245  solvmtx->clustnum, (long)cblk->gcblknum, (long)cblk->ownerid,
246  (long)cblk->lcolidx, (long)colnbr );
247 #endif
248 
249  rc = MPI_Recv( work, colnbr * rhsb->n, PASTIX_MPI_DOUBLE,
250  cblk->ownerid, cblk->gcblknum, solvmtx->solv_comm, &status );
251  assert( rc == MPI_SUCCESS );
252 
253 #if defined (PASTIX_DEBUG_MPI)
254  fprintf( stderr, "[%2d] RHS Fwd: Received cblk %ld from %2d\n",
255  solvmtx->clustnum, (long)cblk->gcblknum, status.MPI_SOURCE );
256 #endif
257 
258  b += cblk->lcolidx;
259  core_dgeadd( PastixNoTrans, colnbr, rhsb->n,
260  1., work, colnbr,
261  1., b, rhsb->ld );
262 
263  (void)rc;
264 #else
265  (void)solvmtx;
266  (void)cblk;
267  (void)work;
268  (void)rhsb;
269 #endif
270 }
271 
272 #if defined( PASTIX_WITH_MPI )
273 /**
274  *******************************************************************************
275  *
276  * @brief Asynchronously sends the rhs associated to cblk->ownerid
277  *
278  *******************************************************************************
279  *
280  * @param[inout] solvmtx
281  * The solver matrix structure.
282  *
283  * @param[in] rhsb
284  * The pointer to the rhs data structure that holds the vectors of the
285  * right hand side.
286  *
287  * @param[in] cblk
288  * The column block that will be sent.
289  *
290  *******************************************************************************/
291 static inline void
292 cpucblk_disend_rhs_fwd( SolverMatrix *solvmtx,
293  pastix_rhs_t rhsb,
294  SolverCblk *cblk )
295 {
296  MPI_Request request;
297  double *b;
298  pastix_int_t colnbr = cblk_colnbr(cblk);
299  pastix_int_t size = colnbr * rhsb->n;
300  pastix_int_t idx = - cblk->bcscnum - 1;
301  int rc;
302 
303  assert( colnbr <= solvmtx->colmax );
304  assert( cblk->cblktype & CBLK_FANIN );
305  assert( rhsb->cblkb[ idx ] != NULL );
306 
307 #if defined(PASTIX_DEBUG_MPI)
308  fprintf( stderr, "[%2d] RHS Fwd: Post Isend cblk %ld to %2d at index %ld of size %ld\n",
309  solvmtx->clustnum, (long)cblk->gcblknum, cblk->ownerid,
310  (long)cblk->lcolidx, (long)(size * sizeof(double) ) );
311 #endif
312 
313  b = (double*)(rhsb->cblkb[ idx ]);
314  assert( b != NULL );
315 
316  rc = MPI_Isend( b, size, PASTIX_MPI_DOUBLE, cblk->ownerid, cblk->gcblknum,
317  solvmtx->solv_comm, &request );
318  assert( rc == MPI_SUCCESS );
319 
320  solverCommMatrixAdd( solvmtx, cblk->ownerid, size * sizeof(double) );
321 
322  /* Register the request to make it progress */
323  pastix_atomic_lock( &(solvmtx->reqlock) );
324 
325  assert( solvmtx->reqidx[ solvmtx->reqnum ] == -1 );
326  assert( solvmtx->reqnum >= 0 );
327  assert( solvmtx->reqnum < solvmtx->reqnbr );
328 
329  solvmtx->reqtab[ solvmtx->reqnum ] = request;
330  solvmtx->reqidx[ solvmtx->reqnum ] = cblk - solvmtx->cblktab;
331  solvmtx->reqnum++;
332 
333  pastix_atomic_unlock( &(solvmtx->reqlock) );
334 
335  (void)rc;
336 }
337 
338 /**
339  *******************************************************************************
340  *
341  * @brief Handle a finished request on a fanin
342  *
343  *******************************************************************************
344  *
345  * @param[in] enums
346  * Enums needed for the solve.
347  *
348  * @param[in] solvmtx
349  * The solver matrix structure.
350  *
351  * @param[in] rhsb
352  * The pointer to the rhs data structure that holds the vectors of the
353  * right hand side.
354  *
355  * @param[inout] cblk
356  * The cblk concerned by the computation.
357  *
358  *******************************************************************************/
359 void
360 cpucblk_drequest_rhs_fwd_handle_send( const args_solve_t *enums,
361  SolverMatrix *solvmtx,
362  pastix_rhs_t rhsb,
363  const SolverCblk *cblk )
364 {
365  pastix_int_t idx = - cblk->bcscnum - 1;
366 
367  assert( cblk->cblktype & CBLK_FANIN );
368  assert( enums->solve_step == PastixSolveForward );
369 
370 #if defined(PASTIX_DEBUG_MPI)
371  {
372  size_t cblksize = cblk_colnbr( cblk ) * rhsb->n * sizeof(double);
373 
374  fprintf( stderr, "[%2d] RHS Fwd: Isend for cblk %ld toward %2d ( %ld Bytes ) (DONE)\n",
375  solvmtx->clustnum, (long)cblk->gcblknum, cblk->ownerid, (long)cblksize );
376  }
377 #endif
378 
379  memFree_null( rhsb->cblkb[ idx ] );
380  (void)solvmtx;
381  (void)enums;
382 }
383 
384 /**
385  *******************************************************************************
386  *
387  * @brief Handle a finished request on a recv cblk.
388  *
389  *******************************************************************************
390  *
391  * @param[in] enums
392  * Enums needed for the solve.
393  *
394  * @param[inout] solvmtx
395  * The solver matrix structure.
396  *
397  * @param[in] rhsb
398  * The pointer to the rhs data structure that holds the vectors of the
399  * right hand side.
400  *
401  * @param[inout] threadid
402  * Id of the thread calling this method.
403  *
404  * @param[inout] status
405  * Statuses of the completed request.
406  *
407  * @param[inout] recvbuf
408  * Reception buffer of the completed request.
409  *
410  *******************************************************************************/
411 static inline void
412 cpucblk_drequest_rhs_fwd_handle_recv( const args_solve_t *enums,
413  SolverMatrix *solvmtx,
414  pastix_rhs_t rhsb,
415  int threadid,
416  const MPI_Status *status,
417  double *recvbuf )
418 {
419  SolverCblk *cblk, *fcbk;
420  pastix_int_t colnbr;
421  int src = status->MPI_SOURCE;
422  int tag = status->MPI_TAG;
423  int nrhs = rhsb->n;
424 
425  assert( ( 0 <= src ) && ( src < solvmtx->clustnbr ) );
426  assert( ( 0 <= tag ) && ( tag < solvmtx->gcblknbr ) );
427 
428  /*
429  * Let's look for the local cblk
430  */
431  fcbk = solvmtx->cblktab + solvmtx->gcbl2loc[ tag ];
432  cblk = fcbk-1;
433 
434  /* Get through source */
435  while( cblk->ownerid != src ) {
436  cblk--;
437  assert( cblk >= solvmtx->cblktab );
438  assert( cblk->gcblknum == tag );
439  assert( cblk->cblktype & CBLK_RECV );
440  }
441  assert( fcbk == (solvmtx->cblktab + cblk->fblokptr->fcblknm) );
442 
443  colnbr = cblk_colnbr( cblk );
444 #if defined(PASTIX_DEBUG_MPI)
445  {
446  int rc;
447  int count = 0;
448  pastix_int_t size = colnbr * nrhs * sizeof(double);
449 
450  rc = MPI_Get_count( status, MPI_CHAR, &count );
451  assert( rc == MPI_SUCCESS );
452  assert( count == size );
453 
454  /* We can't know the sender easily, so we don't print it */
455  fprintf( stderr, "[%2d] RHS Fwd : recv of size %d/%ld for cblk %ld (DONE)\n",
456  solvmtx->clustnum, count, (long)size, (long)cblk->gcblknum );
457  }
458 #endif
459 
460  /* Initialize the cblk with the reception buffer */
461  cblk->threadid = (fcbk->threadid == -1) ? threadid : fcbk->threadid;
462 
463  {
464  double *b = rhsb->b;
465  b += cblk->lcolidx;
466  pastix_cblk_lock( fcbk );
467  core_dgeadd( PastixNoTrans, colnbr, nrhs,
468  1., recvbuf, colnbr,
469  1., b, rhsb->ld );
470  pastix_cblk_unlock( fcbk );
471  }
472 
473  /* Receptions cblks contribute to themselves */
474  cpucblk_drelease_rhs_fwd_deps( enums, solvmtx, rhsb, cblk, fcbk );
475 
476  /* Free the CBLK_RECV */
477  memFree_null( recvbuf );
478 }
479 
480 /**
481  *******************************************************************************
482  *
483  * @brief Handle a finished request.
484  *
485  * If cblktype & CBLK_FANIN : Will deallocate the coeftab
486  * If cblktype & CBLK_RECV : Will add cblk and deallocate the coeftab
487  *
488  *******************************************************************************
489  *
490  * @param[in] enums
491  * Enums needed for the solve.
492  *
493  * @param[inout] solvmtx
494  * The solver matrix structure.
495  *
496  * @param[in] rhsb
497  * The pointer to the rhs data structure that holds the vectors of the
498  * right hand side.
499  *
500  * @param[in] threadid
501  * Id of the thread calling this method.
502  *
503  * @param[in] outcount
504  * Amount of finshed requests.
505  *
506  * @param[in] indexes
507  * Array of completed requests.
508  *
509  * @param[in] statuses
510  * Array of statuses for the completed requests.
511  *
512  *******************************************************************************
513  *
514  * @retval The amount of finished requests updated.
515  *
516  *******************************************************************************/
517 static inline int
518 cpucblk_drequest_rhs_fwd_handle( const args_solve_t *enums,
519  SolverMatrix *solvmtx,
520  pastix_rhs_t rhsb,
521  int threadid,
522  int outcount,
523  const int *indexes,
524  const MPI_Status *statuses )
525 {
526  pastix_int_t i, reqid;
527  int nbrequest = outcount;
528 
529  for ( i = 0; i < outcount; i++ ) {
530  reqid = indexes[i];
531 
532  /*
533  * Handle the reception
534  */
535  if ( solvmtx->reqidx[reqid] == -1 ) {
536  /* We're on a cblk recv, copy datas and restart communications */
537  double *recvbuf;
538  MPI_Status status;
539  int size;
540 
541  memcpy( &status, statuses + i, sizeof(MPI_Status) );
542  MPI_Get_count( &status, PASTIX_MPI_DOUBLE, &size );
543 
544  MALLOC_INTERN( recvbuf, size, double );
545  memcpy( recvbuf, solvmtx->rcoeftab, size * sizeof(double) );
546 
547  solvmtx->recvcnt--;
548 
549  /* Let's restart the communication */
550  assert( solvmtx->recvcnt >= 0 );
551  if ( solvmtx->recvcnt > 0 ) {
552  MPI_Start( solvmtx->reqtab + reqid );
553  nbrequest--;
554  }
555  else {
556  MPI_Request_free( solvmtx->reqtab + reqid );
557  solvmtx->reqtab[reqid] = MPI_REQUEST_NULL;
558  }
559 
560  cpucblk_drequest_rhs_fwd_handle_recv( enums, solvmtx, rhsb,
561  threadid, &status, recvbuf );
562  }
563  /*
564  * Handle the emission
565  */
566  else {
567  SolverCblk *cblk = solvmtx->cblktab + solvmtx->reqidx[ reqid ];
568  assert( cblk->cblktype & CBLK_FANIN );
569 
570  cpucblk_drequest_rhs_fwd_handle_send( enums, solvmtx, rhsb, cblk );
571 
572 #if !defined(NDEBUG)
573  solvmtx->reqidx[ reqid ] = -1;
574 #endif
575  solvmtx->fanincnt--;
576  }
577  }
578 
579  return nbrequest;
580 }
581 
582 /**
583  *******************************************************************************
584  *
585  * @ingroup kernel_fact
586  * @brief Progress communications for one process
587  *
588  * If a communication is completed, it will be treated.
589  * If cblktype & CBLK_FANIN : Will deallocate coeftab
590  * If cblktype & CBLK_RECV : Will add cblk to fcblk
591  *
592  *******************************************************************************
593  *
594  * @param[in] enums
595  * Enums needed for the solve.
596  *
597  * @param[inout] solvmtx
598  * The solver matrix structure.
599  *
600  * @param[in] rhsb
601  * The pointer to the rhs data structure that holds the vectors of the
602  * right hand side.
603  *
604  * @param[in] threadid
605  * Id of the thread calling this method.
606  *
607  *******************************************************************************/
608 void
609 cpucblk_dmpi_rhs_fwd_progress( const args_solve_t *enums,
610  SolverMatrix *solvmtx,
611  pastix_rhs_t rhsb,
612  int threadid )
613 {
614  pthread_t tid = pthread_self();
615  int outcount = 1;
616  int nbrequest, nbfree;
617  int indexes[ solvmtx->reqnbr ];
618  MPI_Status statuses[ solvmtx->reqnbr ];
619 
620  /* Check if someone is already communicating or not */
621  pthread_mutex_lock( &pastix_comm_lock );
622  if ( pastix_comm_tid == (pthread_t)-1 ) {
623  pastix_comm_tid = tid;
624  }
625  pthread_mutex_unlock( &pastix_comm_lock );
626 
627  if ( tid != pastix_comm_tid ) {
628  pastix_yield();
629  return;
630  }
631 
632  /*
633  * Let's register the number of active requests.
634  * We now suppose that the current thread is working on the first nbrequest
635  * active in the reqtab array. Additional requests can be posted during this
636  * progression, but it will be with a larger index. Thus, we do not need to
637  * protect every changes in these requests.
638  * When this is done, the requests arrays is locked to be packed, and the
639  * number of requests is updated for the next round.
640  */
641  pastix_atomic_lock( &(solvmtx->reqlock) );
642  nbrequest = solvmtx->reqnum;
643  pastix_atomic_unlock( &(solvmtx->reqlock) );
644 
645  while( (outcount > 0) && (nbrequest > 0) )
646  {
647  MPI_Testsome( nbrequest, solvmtx->reqtab, &outcount, indexes, statuses );
648  nbfree = 0;
649 
650  /* Handle all the completed requests */
651  if ( outcount > 0 ) {
652  nbfree = cpucblk_drequest_rhs_fwd_handle( enums, solvmtx, rhsb, threadid,
653  outcount, indexes, statuses );
654  }
655 
656  /*
657  * Pack the request arrays, and update the number of active requests by
658  * removing the completed ones
659  */
660  pastix_atomic_lock( &(solvmtx->reqlock) );
661  if ( nbfree > 0 ) {
662  cpucblk_dupdate_reqtab( solvmtx );
663  }
664  nbrequest = solvmtx->reqnum;
665  pastix_atomic_unlock( &(solvmtx->reqlock) );
666  }
667 
668  pastix_comm_tid = (pthread_t)-1;
669  pastix_yield();
670 }
671 #endif /* defined(PASTIX_WITH_MPI) */
672 
673 /**
674  *******************************************************************************
675  *
676  * @brief Wait for incoming dependencies, and return when cblk->ctrbcnt has
677  * reached 0.
678  *
679  *******************************************************************************
680  *
681  * @param[in] rank
682  * The rank of the current thread.
683  *
684  * @param[in] enums
685  * Enums needed for the solve.
686  *
687  * @param[inout] solvmtx
688  * The solver matrix structure.
689  *
690  * @param[inout] cblk
691  * The column block that contribute to fcblk.
692  *
693  * @param[in] rhsb
694  * The pointer to the rhs data structure that holds the vectors of the
695  * right hand side.
696  *
697  *******************************************************************************
698  *
699  * @return 1 if the cblk is a fanin, 0 otherwise
700  *
701  *******************************************************************************/
702 int
704  const args_solve_t *enums,
705  SolverMatrix *solvmtx,
706  SolverCblk *cblk,
707  pastix_rhs_t rhsb )
708 {
709 #if defined(PASTIX_WITH_MPI)
710  if ( cblk->cblktype & CBLK_FANIN ) {
711  /*
712  * We are in the sequential case, we progress on communications and
713  * return if nothing.
714  */
715  //cpucblk_dtestsome( side, solvmtx );
716  return 1;
717  }
718 
719  if ( cblk->cblktype & CBLK_RECV ) {
720  return 1;
721  }
722 
723  /* Make sure we receive every contribution */
724  while( cblk->ctrbcnt > 0 ) {
725  cpucblk_dmpi_rhs_fwd_progress( enums, solvmtx, rhsb, rank );
726  }
727 #else
728  assert( !(cblk->cblktype & (CBLK_FANIN | CBLK_RECV)) );
729  do { pastix_yield(); } while( cblk->ctrbcnt > 0 );
730 #endif
731 
732  (void)rank;
733  (void)enums;
734  (void)solvmtx;
735  (void)rhsb;
736 
737  return 0;
738 }
739 
740 /**
741  *******************************************************************************
742  *
743  * @brief Release the dependencies of the given cblk after an update.
744  *
745  *******************************************************************************
746  *
747  * @param[in] enums
748  * Enums needed for the solve.
749  *
750  * @param[inout] solvmtx
751  * The solver matrix structure.
752  *
753  * @param[in] rhsb
754  * The pointer to the rhs data structure that holds the vectors of the
755  * right hand side.
756  *
757  * @param[in] cblk
758  * The column block that contribute to fcblk.
759  *
760  * @param[inout] fcbk
761  * The facing column block that is updated by cblk.
762  *
763  *******************************************************************************/
764 void
766  SolverMatrix *solvmtx,
767  pastix_rhs_t rhsb,
768  const SolverCblk *cblk,
769  SolverCblk *fcbk )
770 {
771  int32_t ctrbcnt;
772  ctrbcnt = pastix_atomic_dec_32b( &(fcbk->ctrbcnt) );
773  if ( !ctrbcnt ) {
774 #if defined(PASTIX_WITH_MPI)
775  if ( ( fcbk->cblktype & CBLK_FANIN ) &&
776  ( enums->solve_step == PastixSolveForward ) ) {
777  cpucblk_disend_rhs_fwd( solvmtx, rhsb, fcbk );
778  return;
779  }
780 #else
781  (void)enums;
782  (void)rhsb;
783 #endif
784  if ( solvmtx->computeQueue ) {
785  pastix_queue_t *queue = solvmtx->computeQueue[ cblk->threadid ];
786  assert( fcbk->priority != -1 );
787  pqueuePush1( queue, fcbk - solvmtx->cblktab, fcbk->priority );
788  }
789  }
790 }
791 
792 /**
793  *******************************************************************************
794  *
795  * @brief Waitall routine for current cblk request
796  *
797  * It may be possible that some cblk will not be deallocated with the static
798  * scheduler. So a cleanup may be necessary.
799  *
800  *******************************************************************************
801  *
802  * @param[in] enums
803  * Enums needed for the solve.
804  *
805  * @param[in] sched
806  * Define which sched is used
807  * @arg PastixSchedSequential if sequential
808  * @arg PastixSchedStatic if multi-threaded static scheduler
809  * @arg PastixSchedDynamic if multi-threaded dynamic scheduler
810  * No other scheduler is supported.
811  *
812  * @param[inout] solvmtx
813  * The solver matrix structure.
814  *
815  * @param[in] rhsb
816  * The pointer to the rhs data structure that holds the vectors of the
817  * right hand side.
818  *
819  *******************************************************************************/
820 void
822  pastix_int_t sched,
823  SolverMatrix *solvmtx,
824  pastix_rhs_t rhsb )
825 {
826  if ( (sched != PastixSchedSequential) &&
827  (sched != PastixSchedStatic) &&
828  (sched != PastixSchedDynamic) )
829  {
830  return;
831  }
832 #if defined(PASTIX_WITH_MPI)
833  pastix_int_t i;
834  int rc;
835  SolverCblk *cblk;
836  int reqnbr = solvmtx->reqnum;
837  MPI_Status status;
838 
839 #if defined(PASTIX_DEBUG_MPI)
840  fprintf( stderr, "[%2d] Wait for all pending communications\n",
841  solvmtx->clustnum );
842 #endif
843 
844  for ( i=0; i<reqnbr; i++ ) {
845  if ( solvmtx->reqtab[i] == MPI_REQUEST_NULL ) {
846  assert( 0 /* MPI_REQUEST_NULL should have been pushed to the end */ );
847  solvmtx->reqnum--;
848  continue;
849  }
850 
851  /* Make sure that we don't have an already cleaned request in dynamic */
852  assert( solvmtx->reqidx[i] != -1 );
853 
854  rc = MPI_Wait( solvmtx->reqtab + i, &status );
855  assert( rc == MPI_SUCCESS );
856 
857  cblk = solvmtx->cblktab + solvmtx->reqidx[i];
858 
859  /* We should wait only for fanin */
860  assert( cblk->cblktype & CBLK_FANIN );
861 
862  cpucblk_drequest_rhs_fwd_handle_send( enums, solvmtx, rhsb, cblk );
863 
864  solvmtx->reqnum--;
865  }
866  assert( solvmtx->reqnum == 0 );
867  (void)rc;
868 #else
869  (void)enums;
870  (void)solvmtx;
871  (void)rhsb;
872 #endif
873 }
BEGIN_C_DECLS typedef int pastix_int_t
Definition: datatypes.h:51
static void pqueuePush1(pastix_queue_t *q, pastix_int_t elt, double key1)
Push an element with a single key.
Definition: queue.h:64
Queue structure.
Definition: queue.h:38
int core_dgeadd(pastix_trans_t trans, pastix_int_t M, pastix_int_t N, double alpha, const double *A, pastix_int_t LDA, double beta, double *B, pastix_int_t LDB)
Add two matrices together.
Definition: core_dgeadd.c:78
int cpucblk_dincoming_rhs_fwd_deps(int rank, const args_solve_t *enums, SolverMatrix *solvmtx, SolverCblk *cblk, pastix_rhs_t rhsb)
Wait for incoming dependencies, and return when cblk->ctrbcnt has reached 0.
void cpucblk_dsend_rhs_backward(const SolverMatrix *solvmtx, SolverCblk *cblk, pastix_rhs_t rhsb)
Send the rhs associated to a cblk->lcolidx to the remote node.
void cpucblk_drecv_rhs_forward(const SolverMatrix *solvmtx, SolverCblk *cblk, double *work, pastix_rhs_t rhsb)
Receive the rhs associated to a cblk->lcolidx to the remote node.
void cpucblk_drecv_rhs_backward(const SolverMatrix *solvmtx, SolverCblk *cblk, pastix_rhs_t rhsb)
Receive the rhs associated to a cblk->lcolidx to the remote node.
void cpucblk_dsend_rhs_forward(const SolverMatrix *solvmtx, SolverCblk *cblk, pastix_rhs_t rhsb)
Send the rhs associated to a cblk->lcolidx to the remote node.
void cpucblk_drelease_rhs_fwd_deps(const args_solve_t *enums, SolverMatrix *solvmtx, pastix_rhs_t rhsb, const SolverCblk *cblk, SolverCblk *fcbk)
Release the dependencies of the given cblk after an update.
void cpucblk_drequest_rhs_fwd_cleanup(const args_solve_t *enums, pastix_int_t sched, SolverMatrix *solvmtx, pastix_rhs_t rhsb)
Waitall routine for current cblk request.
@ PastixNoTrans
Definition: api.h:445
@ PastixSchedStatic
Definition: api.h:335
@ PastixSchedDynamic
Definition: api.h:338
@ PastixSchedSequential
Definition: api.h:334
void ** cblkb
Definition: pastixdata.h:162
pastix_int_t ld
Definition: pastixdata.h:160
pastix_int_t n
Definition: pastixdata.h:159
Main PaStiX RHS structure.
Definition: pastixdata.h:155
pastix_int_t reqnbr
Definition: solver.h:271
pastix_atomic_lock_t reqlock
Definition: solver.h:273
pastix_int_t priority
Definition: solver.h:183
static pastix_int_t cblk_colnbr(const SolverCblk *cblk)
Compute the number of columns in a column block.
Definition: solver.h:329
pastix_int_t fcblknm
Definition: solver.h:144
MPI_Request * reqtab
Definition: solver.h:269
pastix_int_t recvcnt
Definition: solver.h:217
void * rcoeftab
Definition: solver.h:274
pastix_int_t * gcbl2loc
Definition: solver.h:234
pastix_int_t fanincnt
Definition: solver.h:214
volatile int32_t ctrbcnt
Definition: solver.h:163
SolverBlok * fblokptr
Definition: solver.h:168
volatile int32_t reqnum
Definition: solver.h:272
pastix_int_t gcblknum
Definition: solver.h:174
int threadid
Definition: solver.h:182
pastix_int_t lcolidx
Definition: solver.h:170
pastix_int_t bcscnum
Definition: solver.h:175
SolverCblk *restrict cblktab
Definition: solver.h:228
int8_t cblktype
Definition: solver.h:164
pastix_int_t * reqidx
Definition: solver.h:270
int ownerid
Definition: solver.h:181
Arguments for the solve.
Definition: solver.h:88
Solver column block structure.
Definition: solver.h:161
Solver column block structure.
Definition: solver.h:203