Hitmap 1.3
 All Data Structures Namespaces Files Functions Variables Typedefs Friends Macros Groups Pages
hit_com.c
Go to the documentation of this file.
1 
15 /*
16  * <license>
17  *
18  * Hitmap v1.2
19  *
20  * This software is provided to enhance knowledge and encourage progress in the scientific
21  * community. It should be used only for research and educational purposes. Any reproduction
22  * or use for commercial purpose, public redistribution, in source or binary forms, with or
23  * without modifications, is NOT ALLOWED without the previous authorization of the copyright
24  * holder. The origin of this software must not be misrepresented; you must not claim that you
25  * wrote the original software. If you use this software for any purpose (e.g. publication),
26  * a reference to the software package and the authors must be included.
27  *
28  * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDER AND CONTRIBUTORS "AS IS" AND ANY
29  * EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF
30  * MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL
31  * THE AUTHORS OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
32  * EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
33  * SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
34  * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY,
35  * OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
36  * SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
37  *
38  * Copyright (c) 2007-2015, Trasgo Group, Universidad de Valladolid.
39  * All rights reserved.
40  *
41  * More information on http://trasgo.infor.uva.es/
42  *
43  * </license>
44 */
45 
46 #include <stdio.h>
47 #include <stdlib.h>
48 #include <hit_com.h>
49 #include <unistd.h>
50 #include <hit_funcop.h>
51 #include <hit_sshape.h>
52 #include <hit_cshape.h>
53 #include <hit_bshape.h>
54 #include <hit_error.h>
55 #include <pthread.h>
56 
57 /* Hit NULL CONSTANTS INITIALIZATION */
61 
62 /* Hit DEFAULT REDUCTION OPERATORS */
69 
70 /* Hit EXTRA COMMITED TYPES */
72 
73 
74 MPI_Errhandler mpi_errhandler;
75 
76 
80 static void mpi_error_handler(MPI_Comm *comm, int *err, ... ){
81 
82  HIT_NOT_USED(comm);
83 
84  // 1. Array with space for the MPI Message
85  char array[1024];
86  array[0] = '\0';
87  int lenght;
88 
89  // 2. Check if MPI is initialized to use MPI_Error_string.
90  int initialized;
91  MPI_Initialized(&initialized);
92 
93  if(initialized){
94  MPI_Error_string(*err, array, &lenght);
95  array[lenght] = '\0';
96  }
97 
98 #ifdef __GNUC__
99  // If using gcc we can call the get_gdb_trace function to get
100  // the gdb trace of the current process.
101  printf( "\x1b[31m %s \n \x1b[35m %s \x1b[0m",array,get_gdb_trace());
102 #else
103  printf( "\x1b[31m %s \n",array);
104 #endif
105 
106 }
107 
108 
109 
110 /* Hit MPI INITIALIZATION */
111 void hit_comInit(int *pargc, char **pargv[]) {
112 
113  /* 1. INITIALIZE MPI */
114  int mpi_provided_thread;
115  int ok = MPI_Init_thread( pargc, pargv, MPI_THREAD_MULTIPLE, &mpi_provided_thread );
116  hit_mpiTestError( ok, "Initializing HitCom for MPI comm. library" );
117 
118  // Set our error handler
119  MPI_Comm_create_errhandler(mpi_error_handler, &mpi_errhandler);
120  MPI_Comm_set_errhandler(MPI_COMM_WORLD, mpi_errhandler);
121 
122  /* 2. GET BASIC GROUP PARAMETERS */
123  hit_malloc( HIT_TOPOLOGY_INFO, HitPTopology, 1 );
124  // @arturo Ago 2015: New allocP interface
125  // sizeof( HitPTopology), HitPTopology*);
126  MPI_Comm_rank(MPI_COMM_WORLD, &(HIT_TOPOLOGY_INFO->selfRank) );
127  MPI_Comm_size(MPI_COMM_WORLD, &(HIT_TOPOLOGY_INFO->numProcs) );
128  MPI_Comm_dup ( MPI_COMM_WORLD, &(HIT_TOPOLOGY_INFO->comm) );
129  HIT_TOPOLOGY_INFO->numUses = 1;
130  // @arturo 2015/01/03,
131  // New functions to redistribute need proper initialization of fields "next" and new "global"
132  HIT_TOPOLOGY_INFO->next = NULL;
134  // @arturo 2015/01/07
135  HIT_TOPOLOGY_INFO->antiCommGroup = MPI_GROUP_EMPTY;
136 
137  /* 3. SET SOME NULL VALUES TO MPI EQUIVALENTS */
138  HIT_RANK_NULL = MPI_PROC_NULL;
139  int i;
140  for (i=0; i<HIT_MAXDIMS; i++) HIT_RANKS_NULL.rank[i] = MPI_PROC_NULL;
141 
142  /* 4. DECLARE THE BASIC REDUCTION OPERATIONS ON int AND double */
149 
150  /* 5. COMMIT EXTRA DATATYPES */
151  ok = MPI_Type_contiguous( HIT_MAXDIMS * HIT_SIG_NUM_FIELDS + 1, MPI_INT, &HIT_SHAPE_SIG );
152  hit_mpiTestError(ok, "Creating shape datatype");
153  ok = MPI_Type_commit( &HIT_SHAPE_SIG );
154  hit_mpiTestError(ok, "Commiting shape datatype");
155 }
156 
157 
158 /* Hit MPI FINALIZATION */
160 
161  /* 1. FREE BASIC TOPOLOGY INFORMATION */
162  //MPI_Comm_free( (MPI_Comm *) HIT_TOPOLOGY_INFO.lowLevel);
163  //free(HIT_TOPOLOGY_INFO.lowLevel);
164  //free(HIT_TOPOLOGY_INFO.lowLevel);
165  if ( HIT_TOPOLOGY_INFO->numUses > 1 )
166  hit_warnInternal( __FUNCTION__, "Not all topologies or layouts have been deallocted", "", __FILE__, __LINE__);
169 
170  /* 2. FREE BASIC REDUCTION OPERATIONS ON int AND double */
177 
178  /* 3. FREE COMMITED EXTRA DATATYPES */
179  int ok = MPI_Type_free( &HIT_SHAPE_SIG );
180  hit_mpiTestError(ok, "Freeing shape datatype");
181 
182  /* 3. FINALIZE MPI */
183  MPI_Finalize();
184  pthread_exit(NULL);
185 }
186 
187 
188 /* Hit COMPUTE MPI DERIVED DATA TYPE FOR A SINGLE TILE */
189 HitType hit_comType(const void *varP, HitType baseType) {
190  HitType newType, prevType;
191  HitType aux[HIT_MAXDIMS];
192  int dim,numTypes, contiguous;
193  MPI_Aint baseExtent;
194  HitTile var = *(const HitTile *)varP;
195 
196 #ifdef BARELY_FASTER
197  /* 1. ORIGINAL VARIABLES OR SELECTIONS OF THE WHOLE SHAPE */
198  if (var.acumCard == var.origAcumCard[0]) {
199 
200  /* 1.2. CONTIGUOS TYPE */
201  if (MPI_SUCCESS != MPI_Type_contiguous(var.acumCard, baseType, &newType))
202  hit_error("MPI_Type_contiguous", __FILE__, __LINE__);
203  if (MPI_SUCCESS != MPI_Type_commit(&newType))
204  hit_error("MPI_Type_commit", __FILE__, __LINE__);
205 
206  /* 1.3. RETURN */
207  return newType;
208  }
209 #endif
210 
211  /* 2. GET MEMORY FOR THE TYPES ARRAY */
212  numTypes=hit_tileDims(var);
213 
214  /* 3. GET EXTENT OF THE BASE TYPE FOR THE STRIDES */
215  MPI_Type_extent(baseType, &baseExtent);
216 
217  /* 4. START LOOP OF DIMENSIONS BACKWARDS (ROW ORIENTED) */
218  prevType = baseType;
219  contiguous = 1;
220  int partialAcumCard = 1;
221 
222  for (dim=numTypes-1; dim >=0; dim-- ) {
223  /* 4.1. WHILE CONTIGUOUS AND THIS DIM HAS NO STRIDE, KEEP ON CONTIGUOS DATA TYPES */
224  if (contiguous && hit_tileDimStride(var,dim) == 1 ) {
225 #ifdef DEBUG
226 printf("Contiguous type for dim %d, with %d total elements\n", dim, partialAcumCard); fflush(stdout);
227 #endif
228  /* 4.1.1. COMPUTE TYPE */
229  if (MPI_SUCCESS != MPI_Type_contiguous( hit_tileDimCard(var,dim), prevType,
230  &(aux[dim])))
231  hit_error("MPI_Type_contiguous", __FILE__, __LINE__);
232 
233  /* 4.1.2. A NON-FULL DIMENSION FORCES TO USE HVECTOR TYPES AFTERWARDS */
234  partialAcumCard = partialAcumCard * hit_tileDimCard(var, dim);
235  if ( partialAcumCard != var.origAcumCard[dim]) contiguous = 0;
236  }
237  /* 4.2. AFTER A NON-CONTIGUOS DIMENSION, USE GENERIC HVECTOR TYPES */
238  else {
239  contiguous = 0;
240 #ifdef DEBUG
241 printf("Hvector type for dim %d, Card: %d, stride: %d\n", dim, var.card[dim], baseExtent*var.origAcumCard[dim+1]* hit_tileDimStride(var,dim)
242 //.shape.sig[dim].stride
243 ); fflush(stdout);
244 #endif
245  if (MPI_SUCCESS != MPI_Type_create_hvector(var.card[dim], 1,
246  baseExtent*var.origAcumCard[dim+1]* hit_tileDimSig(var,dim).stride,
247  prevType, &(aux[dim])) )
248  hit_error("MPI_Type_create_hvector", __FILE__, __LINE__);
249  }
250  prevType = aux[dim];
251  }
252 
253 #ifdef DEBUG
254  printf("freeing type resources\n"); fflush(stdout);
255 #endif
256  newType=prevType;
257 
258  for (dim=numTypes-1; dim >0; dim-- ) {
259  hit_comFreeType(aux[dim]);
260  }
261 
262 #ifdef DEBUG
263  printf("commited hit type %d\n",newType);
264 #endif
265 
266  /* 6. RETURN THE TYPES STRUCTURE */
267  return newType;
268 }
269 
278 inline void * hit_comSearchData (const void * varP) {
279  const HitTile * var = (const HitTile *)varP;
280  while ( var->hierDepth > 0 ) {
281  var=(HitTile *)(var->data);
282  }
283  return var->data;
284 }
285 
286 /*
287  * CREATE A DERIVED DATA TYPE FOR A TILE
288  * RECURSIVELY FOR MULTI-LEVEL TILES
289  */
290 HitType hit_comTypeRec(const void *varP, HitType baseType) {
291  HitTile var = *(const HitTile *)varP;
292  HitType *types;
293 
294  MPI_Aint *addresses;
295  int i, j, offset, *sizes, ind[HIT_MAXDIMS], max[HIT_MAXDIMS];
296  HitTile * auxtile;
297  HitType newType;
298  HitType * auxtypes;
299 
300  if ( var.memStatus == HIT_MS_NULL || hit_tileDims(var)==0 ){
301 #ifdef DEBUG
302  printf("hit_comTypeRec: NULL type, Tile is NULL or dims == 0\n");
303 #endif
304  return HIT_TYPE_NULL;
305  }
306 
307  if( var.hierDepth > 0 ) {
308  hit_malloc(sizes, int, var.acumCard );
309  hit_malloc(addresses, MPI_Aint, var.acumCard );
310  hit_malloc(types, HitType, var.acumCard );
311  hit_malloc(auxtypes, HitType, var.acumCard );
312  // @arturo Ago 2015: New allocP interface
313  /*
314  hit_malloc(sizes, int, (size_t)var.acumCard*sizeof(int),int*);
315  hit_malloc(addresses,(size_t)var.acumCard*sizeof(MPI_Aint),MPI_Aint*);
316  hit_malloc(types,(size_t)var.acumCard*sizeof(HitType),HitType*);
317  hit_malloc(auxtypes,(size_t)var.acumCard*sizeof(HitType),HitType*);
318  */
319  //Initialize max with cardinalities
320  memcpy(max,var.card,(size_t) hit_tileDims(var) *sizeof(int));
321  //initialize ind with 0
322  memset(ind,0,(size_t)hit_tileDims(var)*sizeof(int));
323 #ifdef DEBUG
324 printf("Multilevel Struct type for acumCard: %d\n", var.acumCard ); fflush(stdout);
325 #endif
326  for(i=0;i<var.acumCard;i++) {
327  offset=0;
328  for(j=0;j<hit_tileDims(var);j++) {
329  offset+=ind[j]*var.qstride[j]*var.origAcumCard[j+1];
330  }
331  auxtile=(HitTile *)( (char *)var.data+offset*(int)var.baseExtent );
332  auxtypes[i]=hit_comTypeRec(auxtile, baseType);
333  types[i]=auxtypes[i];
334  void * data = hit_comSearchData(auxtile);
335  MPI_Get_address(data,&addresses[i]);
336  sizes[i]=1;
337  j=0;
338  do {
339  ind[j]=(ind[j]+1)%max[j];
340  j++;
341  } while(j<hit_tileDims(var) && ind[j-1]==0 && max[j]!=0);
342  }
343  for(i=var.acumCard-1;i>=0;i--) addresses[i]-=addresses[0];
344  MPI_Type_create_struct(var.acumCard,sizes,addresses,types,&newType);
345  MPI_Type_commit(&newType);
346  for(i=0;i<var.acumCard;i++) hit_comFreeType(auxtypes[i]);
347  free(auxtypes);
348  free(types);
349  free(sizes);
350  free(addresses);
351  return newType;
352  } else {
353  newType = hit_comType(varP,baseType);
354  /* COMMIT TYPE */
355  if (MPI_SUCCESS != MPI_Type_commit(&newType))
356  hit_error("Commit last type", __FILE__, __LINE__);
357 
358 #ifdef DEBUG
359  printf("commited hit type %d\n",newType);
360 #endif
361  return newType;
362  }
363 }
364 
365 
366 /* Create the low level communicator to allow collective dim communications */
367 void hit_comAllowDim(HitLayout * lay, int dim){
368 
369  // @arturo Mar 2013
370  HitPTopology *ptopo = NULL;
371 
372  if(lay->active){
373  // Get the active ranks for the processor
374  HitRanks activeRanks = hit_laySelfRanks( *lay );
375 
376  // Create the color to split the low level communicator
377  int colour=0, i, acum = 1;
378  for (i=0;i<hit_layNumDims(*lay);i++) {
379  if(i!=dim) colour+= acum * activeRanks.rank[i];
380  acum *= lay->numActives[i];
381  }
382 
383  // Split the low level communicator
384  ptopo = hit_ptopSplit( lay->pTopology[0], colour );
385 
386  }
387 
388  lay->pTopology[ dim+1 ] = ptopo;
389 }
390 
391 
392 
394  HitRanks sendTo,
395  const void *tileP,
396  HitShape selection,
397  int selectMode,
398  HitRanks receiveFrom,
399  HitType baseType,
400  int tag){
401 
402  /* 1. AVOID COMMUNICATION ON INACTIVE PROCS IN THE LAYOUT */
403  if ( ! lay.active ) return HIT_COM_NULL;
404 
405  /* 2. NO COMMUNICATIONS TO/FROM NULL */
406  // @arturo: CORRECTED BUG: Let shifts across boundaries be transformed in communications
407  // from/to NULL ranks, that means from/to MPI_NULL_PROC identifier
408  // Thus, a boundary processor may send but not receive, or viceversa.
409  //if ( hit_ranksCmp(sendTo,HIT_RANKS_NULL) || hit_ranksCmp(receiveFrom,HIT_RANKS_NULL) )
410  // return HIT_COM_NULL;
411 
412  /* 3. GET TOPOLOGY COMMUNICATOR */
413  //MPI_Comm comm = *((MPI_Comm *)lay.topo.pTopology.lowLevel);
414  MPI_Comm comm = lay.topo.pTopology->comm;
415 
416  /* 4. DECLARE NEW COMMUNICATION ISSUE */
417  HitCom newCom = HIT_COM_NULL;
418 
419  /* 5. IF TILE HAS NO MEMORY OR NO DIMENSIONS RETURN NULL COMMUNICATION */
420  HitTile Select;
421  const HitTile *tile = (const HitTile *)tileP;
422  if ( tile->memStatus == HIT_MS_NULL || tile->memStatus == HIT_MS_NOMEM || hit_tileDims(*tile)==0 ) return HIT_COM_NULL;
423 
424  /* 6. SELECTION */
425  if (selectMode == HIT_COM_TILECOORDS) hit_tileSelect( &Select, tile, selection );
426  else if (selectMode == HIT_COM_ARRAYCOORDS) hit_tileSelectArrayCoords( &Select, tile, selection);
427  else return HIT_COM_NULL;
428 
429  /* 7. TYPES */
430  newCom.typeSend = hit_comTypeRec( &Select, baseType);
431  newCom.dataSend = hit_comSearchData(&Select);
432 
433  /* 8. FILLING UP OTHER FIELDS */
434  // @arturo 2015/01/22
435  //newCom.myself = lay.topo.linearRank;
436  newCom.myself = hit_topSelfRankInternal( lay.topo );
437  newCom.comm = comm;
438  newCom.tag = tag;
440 
441  newCom.sendTo = hit_topRankInternal(lay.topo,sendTo);
442  newCom.recvFrom = hit_topRankInternal(lay.topo,receiveFrom);
443 
444  /* 9. RETURN */
445  return newCom;
446 }
447 
448 
449 
451  HitRanks sendTo,
452  const void *tilePSend,
453  HitShape selectionSend,
454  int selectSendMode,
455  HitRanks receiveFrom,
456  const void *tilePRecv,
457  HitShape selectionRecv,
458  int selectRecvMode,
459  HitType baseType,
460  int tag){
461 
462  /* @arturo Feb 2013: BUG CORRECTED
463  * Processors that are inactive in the layout may send/receive data
464  * In hit_patternRedistribute there are two layouts involved, sometimes
465  * the local process is not active in one layout, but it is in the other.
466  */
467  /* 1. AVOID COMMUNICATION ON INACTIVE PROCS IN THE TOPOLOGY */
468  //if ( ! lay.active ) return HIT_COM_NULL;
469  if ( ! lay.topo.active ) return HIT_COM_NULL;
470 
471  /* 1. GET TOPOLOGY COMMUNICATOR */
472  //MPI_Comm comm = *((MPI_Comm *)lay.topo.pTopology.lowLevel);
473  MPI_Comm comm = lay.topo.pTopology->comm;
474 
475  /* 2. DECLARE NEW COMMUNICATION ISSUE */
476  HitCom newCom = HIT_COM_NULL;
477 
478  /* 3. START FILLING UP ITS FIELDS */
479  // @arturo 2015/01/22
480  //newCom.myself = lay.topo.linearRank;
481  newCom.myself = hit_topSelfRankInternal( lay.topo );
482 
483  /* 4.1. IGNORE COMMUNICATIONS TO NULL */
484  /* TO THE SAME PROCESSOR ARE NEEDED IN CASE THE SEND/RECV BUFFERS ARE DIFFERENT */
485  if ( hit_ranksCmp(sendTo,HIT_RANKS_NULL) ) {
486 #ifdef DEBUG
487  printf("%d no sending\n", hit_Rank); fflush(stdout);
488 #endif
489  newCom.sendTo = MPI_PROC_NULL;
490  newCom.typeSend = HIT_TYPE_NULL;
491  newCom.dataSend = NULL;
492  }
493  /* 4.2. COMMUNICATIONS TO MYSELF OR ANOTHER PROCCESSOR */
494  else {
495  newCom.sendTo = hit_topRankInternal(lay.topo,sendTo);
496 #ifdef DEBUG
497  printf("[%d] Create Send Selection -- SendTo: %d\n", hit_Rank, newCom.sendTo); fflush(stdout);
498 #endif
499  const HitTile *tileSend = (const HitTile *)tilePSend;
500  HitTile selectSend;
501  /* IF TILE HAS NO MEMORY OR NO DIMENSIONS RETURN NULL COMMUNICATION */
502  if ( tileSend->memStatus == HIT_MS_NULL || tileSend->memStatus == HIT_MS_NOMEM || hit_tileDims(*tileSend)==0 ) return HIT_COM_NULL;
503 
504  /* @arturo: EMPTY SELECTION LEADS TO NULL TYPE, NO COMMUNICATION */
505  if((selectSendMode == HIT_COM_TILECOORDS && !hit_tileCheckBoundaryTileCoords(tileSend, selectionSend)) ||
506  (selectSendMode == HIT_COM_ARRAYCOORDS && !hit_tileCheckBoundaryArrayCoords(tileSend, selectionSend))
507  ){
508 #ifdef DEBUG
509  printf("send select out of boundary\n"); fflush(stdout);
510 #endif
511  return HIT_COM_NULL;
512  }
513 
514 
515  if (selectSendMode == HIT_COM_TILECOORDS) hit_tileSelect( &selectSend, tileSend, selectionSend );
516  else if (selectSendMode == HIT_COM_ARRAYCOORDS) hit_tileSelectArrayCoords( &selectSend, tileSend, selectionSend);
517  else return HIT_COM_NULL;
518 
519  //printf("Tile %p(%d), Select Tile %p(%d)\n",tileSend,tileSend->memStatus,&selectSend,selectSend.memStatus);
520 
521 
522  newCom.typeSend = hit_comTypeRec( &selectSend, baseType);
523  newCom.dataSend = hit_comSearchData(&selectSend);
524  }
525 
526  /* 4.4. IGNORE COMMUNICATIONS FROM NULL */
527  /* TO THE SAME PROCESSOR ARE NEEDED IN CASE THE SEND/RECV BUFFERS ARE DIFFERENT */
528  if ( hit_ranksCmp(receiveFrom,HIT_RANKS_NULL) ) {
529 #ifdef DEBUG
530  printf("%d no reception\n", hit_Rank); fflush(stdout);
531 #endif
532  newCom.recvFrom = MPI_PROC_NULL;
533  newCom.typeRecv = HIT_TYPE_NULL;
534  newCom.dataRecv = NULL;
535  }
536  else {
537  newCom.recvFrom = hit_topRankInternal(lay.topo,receiveFrom);
538 #ifdef DEBUG
539  printf("[%d] Create Recv Selection -- RecvFrom: %d\n", hit_Rank, newCom.recvFrom); fflush(stdout);
540 #endif
541  const HitTile *tileRecv = (const HitTile *)tilePRecv;
542  HitTile selectRecv;
543  /* IF TILE HAS NO MEMORY OR NO DIMENSIONS RETURN NULL COMMUNICATION */
544  if ( tileRecv->memStatus == HIT_MS_NULL || tileRecv->memStatus == HIT_MS_NOMEM || hit_tileDims(*tileRecv)==0 ) return HIT_COM_NULL;
545 
546  /* @arturo: EMPTY SELECTION LEADS TO NULL TYPE, NO COMMUNICATION */
547 
548  if((selectSendMode == HIT_COM_TILECOORDS && !hit_tileCheckBoundaryTileCoords(tileRecv, selectionSend)) ||
549  (selectSendMode == HIT_COM_ARRAYCOORDS && !hit_tileCheckBoundaryArrayCoords(tileRecv, selectionSend))
550  ){
551 #ifdef DEBUG
552  printf("recv select out of boundary\n"); fflush(stdout);
553 #endif
554  return HIT_COM_NULL;
555  }
556 
557  if (selectRecvMode == HIT_COM_TILECOORDS) hit_tileSelect( &selectRecv, tileRecv, selectionRecv );
558  else if (selectRecvMode == HIT_COM_ARRAYCOORDS) hit_tileSelectArrayCoords( &selectRecv, tileRecv, selectionRecv);
559  else return HIT_COM_NULL;
560 
561  newCom.typeRecv = hit_comTypeRec( &selectRecv, baseType);
562  newCom.dataRecv = hit_comSearchData(&selectRecv);
563  }
564 #ifdef DEBUG
565  printf("[%d] communication created\n", hit_Rank); fflush(stdout);
566 #endif
567 
568  /* 5. FILL UP THE OTHER FIELDS */
569  newCom.comm = comm;
570  newCom.tag = tag;
571  newCom.requestSend = MPI_REQUEST_NULL;
572  newCom.commType = HIT_SENDRECV;
573  /*
574  newCom.statusSend = ???;
575  newCom.statusRecv = ???;
576  */
577 
578  /* 6. RETURN */
579  return newCom;
580 }
581 
582 
583 
584 /* Hit COMMIT REDUCE */
586  HitRanks root,
587  const void * tilePSend,
588  HitShape selectionSend,
589  int selectSendMode,
590  const void * tilePRecv,
591  HitShape selectionRecv,
592  int selectRecvMode,
593  HitType baseType,
594  HitOp operation){
595 
596 
597  /* 2. DECLARE NEW COMMUNICATION ISSUE */
598  HitCom newCom = HIT_COM_NULL;
599 
600  /* 3. GET TOPOLOGY COMMUNICATOR */
601  //MPI_Comm comm = *((MPI_Comm *)lay.pTopology[0].lowLevel);
602  MPI_Comm comm = lay.pTopology[0]->comm;
603  newCom.comm = comm;
604 
605  /* 4. START FILLING UP ITS FIELDS */
606  // @arturo 2015/01/22
607  //newCom.myself = lay.topo.linearRank;
608  newCom.myself = hit_topSelfRankInternal( lay.topo );
609 
610  /* 5. COMMUNICATION */
611  /* 5.1 SEND */
612  HitRanks rootActive = hit_layToActiveRanks(lay,root);
613  newCom.sendTo = hit_layActiveRanksId(lay,rootActive);
614 
615  const HitTile *tileSend = (const HitTile *)tilePSend;
616  HitTile selectSend;
617 
618  // @arturo, Feb 2013
619  // HIT_TILE_NULL in the send tile parameter indicates MPI_IN_PLACE
620  if ( tileSend->memStatus == HIT_MS_NULL ) {
621  newCom.dataSend = MPI_IN_PLACE;
622  }
623  else {
624  if (selectSendMode == HIT_COM_TILECOORDS) hit_tileSelect( &selectSend, tileSend, selectionSend );
625  else if (selectSendMode == HIT_COM_ARRAYCOORDS) hit_tileSelectArrayCoords( &selectSend, tileSend, selectionSend);
626  else return HIT_COM_NULL;
627 
628  newCom.dataSend = hit_comSearchData(&selectSend);
629  }
630 
631  /* 5.2 RECV */
632  const HitTile *tileRecv = (const HitTile *)tilePRecv;
633  HitTile selectRecv;
634 
635  /* if tile has no memory or no dimensions return null communication */
636  if ( tileRecv->memStatus == HIT_MS_NULL || tileRecv->memStatus == HIT_MS_NOMEM || hit_tileDims(*tileRecv)==0 ) return HIT_COM_NULL;
637 
638  if (selectRecvMode == HIT_COM_TILECOORDS) hit_tileSelect( &selectRecv, tileRecv, selectionRecv );
639  else if (selectRecvMode == HIT_COM_ARRAYCOORDS) hit_tileSelectArrayCoords( &selectRecv, tileRecv, selectionRecv);
640  else return HIT_COM_NULL;
641 
642  newCom.dataRecv = hit_comSearchData(&selectRecv);
643 
644  /* 6. COMMUNICATION TYPE */
645  newCom.commType = HIT_REDUCE;
646 
647  /* 7. REDUCTION OPERATION */
648  newCom.operation = operation;
649 
650  /* 8. TYPE */
651  // @arturo, Feb 2013, Due to the new MPI_IN_PLACE option, the selection of send may not exist.
652  // Use recv selection to compute the type
653  //newCom.typeSend = hit_comTypeRec( &selectSend, baseType);
654  newCom.typeSend = hit_comTypeRec( &selectRecv, baseType);
655 
656  /* 9. RETURN */
657  return newCom;
658 }
659 
661  int dim,
662  HitRanks root,
663  const void * tilePSend,
664  HitShape selectionSend,
665  int selectSendMode,
666 
667  const void * tilePRecv,
668  HitShape selectionRecv,
669  int selectRecvMode,
670  HitType baseType,
671  HitOp operation){
672 
673 
674  /* 2. DECLARE NEW COMMUNICATION ISSUE */
675  HitCom newCom = HIT_COM_NULL;
676 
677  /* 3. GET TOPOLOGY COMMUNICATOR */
678  // @arturo Mar 2013
679  //if (lay.pTopology[dim+1].lowLevel==NULL) return HIT_COM_NULL;
680  //MPI_Comm comm = *((MPI_Comm *)lay.pTopology[dim+1].lowLevel);
681  if ( lay.pTopology[dim+1] == NULL ) return HIT_COM_NULL;
682  MPI_Comm comm = lay.pTopology[ dim+1 ]->comm;
683  newCom.comm = comm;
684 
685  /* 4. START FILLING UP ITS FIELDS */
686  // @arturo 2015/01/22
687  //newCom.myself = lay.topo.linearRank;
688  newCom.myself = hit_topSelfRankInternal( lay.topo );
689 
690  /* 5. COMMUNICATION */
691  /* 5.1 SEND */
692  newCom.sendTo = hit_topRankInternal(lay.topo,root);
693  const HitTile *tileSend = (const HitTile *)tilePSend;
694  HitTile selectSend;
695  /* IF TILE HAS NO MEMORY OR NO DIMENSIONS RETURN NULL COMMUNICATION */
696  if ( tileSend->memStatus == HIT_MS_NULL || tileSend->memStatus == HIT_MS_NOMEM || hit_tileDims(*tileSend)==0 ) return HIT_COM_NULL;
697 
698  if (selectSendMode == HIT_COM_TILECOORDS) hit_tileSelect( &selectSend, tileSend, selectionSend );
699  else if (selectSendMode == HIT_COM_ARRAYCOORDS) hit_tileSelectArrayCoords( &selectSend, tileSend, selectionSend);
700  else return HIT_COM_NULL;
701 
702  newCom.dataSend = hit_comSearchData(&selectSend);
703 
704  /* 5.2 RECV */
705  const HitTile *tileRecv = (const HitTile *)tilePRecv;
706  HitTile selectRecv;
707  /* IF TILE HAS NO MEMORY OR NO DIMENSIONS RETURN NULL COMMUNICATION */
708  if ( tileRecv->memStatus == HIT_MS_NULL || tileRecv->memStatus == HIT_MS_NOMEM || hit_tileDims(*tileRecv)==0 ) return HIT_COM_NULL;
709 
710  if (selectRecvMode == HIT_COM_TILECOORDS) hit_tileSelect( &selectRecv, tileRecv, selectionRecv );
711  else if (selectRecvMode == HIT_COM_ARRAYCOORDS) hit_tileSelectArrayCoords( &selectRecv, tileRecv, selectionRecv);
712  else return HIT_COM_NULL;
713 
714  newCom.dataRecv = hit_comSearchData(&selectRecv);
715 
716  /* 6. COMMUNICATION TYPE */
717  newCom.commType = HIT_REDUCE;
718 
719  /* 7. REDUCTION OPERATION */
720  newCom.operation = operation;
721 
722  /* 8. TYPE */
723  newCom.typeSend = hit_comTypeRec( &selectSend, baseType);
724 
725  /* 9. RETURN */
726  return newCom;
727 }
728 
729 
730 /* Hit COMMIT BROADCAST */
732  HitRanks root,
733  const void * tileP,
734  HitShape selection,
735  int selectMode,
736  HitType baseType){
737 
738 
739  /* 2. DECLARE NEW COMMUNICATION ISSUE */
740  HitCom newCom = HIT_COM_NULL;
741 
742  /* 3. GET TOPOLOGY COMMUNICATOR */
743  //MPI_Comm comm = *((MPI_Comm *)lay.pTopology[0].lowLevel);
744  MPI_Comm comm = lay.pTopology[0]->comm;
745  newCom.comm = comm;
746 
747  /* 4. START FILLING UP ITS FIELDS */
748  // @arturo 2015/01/22
749  //newCom.myself = lay.topo.linearRank;
750  newCom.myself = hit_topSelfRankInternal( lay.topo );
751 
752  /* 5. COMMUNICATION */
753  /* 5.1 SEND */
754  HitRanks rootActive = hit_layToActiveRanks(lay,root);
755  newCom.sendTo = hit_layActiveRanksId(lay,rootActive);
756 
757  const HitTile *tile = (const HitTile *)tileP;
758  HitTile Select;
759  /* IF TILE HAS NO MEMORY OR NO DIMENSIONS RETURN NULL COMMUNICATION */
760  if ( tile->memStatus == HIT_MS_NULL || tile->memStatus == HIT_MS_NOMEM || hit_tileDims(*tile)==0 ) return HIT_COM_NULL;
761 
762  if (selectMode == HIT_COM_TILECOORDS) hit_tileSelect( &Select, tile, selection );
763  else if (selectMode == HIT_COM_ARRAYCOORDS) hit_tileSelectArrayCoords( &Select, tile, selection);
764  else return HIT_COM_NULL;
765 
766  newCom.dataSend = hit_comSearchData(&Select);
767 
768  /* 6. COMMUNICATION TYPE */
769  newCom.commType = HIT_BROADCAST;
770 
771  /* 7. TYPE */
772  newCom.typeSend = hit_comTypeRec( &Select, baseType);
773 
774  /* 8. RETURN */
775  return newCom;
776 }
777 
778 
780  int dim,
781  int root,
782  const void * tileP,
783  HitShape selection,
784  int selectMode,
785  HitType baseType){
786 
787  /* 1. IF TILE HAS NO MEMORY OR NO DIMENSIONS RETURN NULL COMMUNICATION */
788  const HitTile *tile = (const HitTile *)tileP;
789  if ( tile->memStatus == HIT_MS_NULL || tile->memStatus == HIT_MS_NOMEM || hit_tileDims(*tile)==0 ) return HIT_COM_NULL;
790 
791  /* 2. DECLARE NEW COMMUNICATION ISSUE */
792  HitCom newCom = HIT_COM_NULL;
793 
794  /* 3. GET TOPOLOGY COMMUNICATOR */
795  // @arturo Mar 2013
796  //if (lay.pTopology[dim+1].lowLevel==NULL) return HIT_COM_NULL;
797  //MPI_Comm comm = *((MPI_Comm *)lay.pTopology[dim+1].lowLevel);
798  if ( lay.pTopology[dim+1] == NULL ) return HIT_COM_NULL;
799  MPI_Comm comm = lay.pTopology[ dim+1 ]->comm;
800  newCom.comm = comm;
801 
802  /* 4. START FILLING UP ITS FIELDS */
803  // @arturo 2015/01/22
804  //newCom.myself = lay.topo.linearRank;
805  newCom.myself = hit_topSelfRankInternal( lay.topo );
806 
807  /* 5. COMMUNICATION */
808  /* 5.1 SEND */
809  HitRanks rootRanks = hit_laySelfRanks( lay );
810  if ( root != HIT_COM_MYSELF ) rootRanks.rank[ dim ] = root;
811  HitRanks rootActive = hit_layToActiveRanks( lay, rootRanks );
812  newCom.sendTo = rootActive.rank[dim];
813  HitTile Select;
814  if (selectMode == HIT_COM_TILECOORDS) hit_tileSelect( &Select, tile, selection );
815  else if (selectMode == HIT_COM_ARRAYCOORDS) hit_tileSelectArrayCoords( &Select, tile, selection);
816  else {
817  fprintf(stderr,"Hit Warning: hit_comBroadcastDimSelect: Unknown constant for coordinate type\n");
818  return HIT_COM_NULL;
819  }
820 
821  newCom.dataSend = hit_comSearchData(&Select);
822 
823  /* 6. COMMUNICATION TYPE */
824  newCom.commType = HIT_BROADCAST;
825 
826  /* 7. DATA TYPE */
827  newCom.typeSend = hit_comTypeRec( &Select, baseType);
828 
829  /* 8. RETURN */
830  return newCom;
831 }
832 
833 
835  const void * tilePSend,
836  HitShape selectionSend,
837  int selectSendMode,
838  const void * tilePRecv,
839  HitShape selectionRecv,
840  int selectRecvMode,
841  HitType baseType,
842  int count){
843 
844 
845  /* 2. DECLARE NEW COMMUNICATION ISSUE */
846  HitCom newCom = HIT_COM_NULL;
847 
848  /* 3. GET TOPOLOGY COMMUNICATOR */
849  //MPI_Comm comm = *((MPI_Comm *)lay.pTopology[0].lowLevel);
850  MPI_Comm comm = lay.pTopology[0]->comm;
851  newCom.comm = comm;
852 
853  /* 4. START FILLING UP ITS FIELDS */
854  // @arturo 2015/01/22
855  //newCom.myself = lay.topo.linearRank;
856  newCom.myself = hit_topSelfRankInternal( lay.topo );
857 
858  /* 5. COMMUNICATION */
859  /* 5.1. SEND */
860  const HitTile *tileSend = (const HitTile *)tilePSend;
861  HitTile selectSend;
862  /* IF TILE HAS NO MEMORY OR NO DIMENSIONS RETURN NULL COMMUNICATION */
863  if ( tileSend->memStatus == HIT_MS_NULL || tileSend->memStatus == HIT_MS_NOMEM || hit_tileDims(*tileSend)==0 ) return HIT_COM_NULL;
864 
865  if (selectSendMode == HIT_COM_TILECOORDS) hit_tileSelect( &selectSend, tileSend, selectionSend );
866  else if (selectSendMode == HIT_COM_ARRAYCOORDS) hit_tileSelectArrayCoords( &selectSend, tileSend, selectionSend);
867  else return HIT_COM_NULL;
868 
869  newCom.dataSend = selectSend.data;
870 
871  /* 5.2. RECV */
872  const HitTile *tileRecv = (const HitTile *)tilePRecv;
873  HitTile selectRecv;
874  /* IF TILE HAS NO MEMORY OR NO DIMENSIONS RETURN NULL COMMUNICATION */
875  if ( tileRecv->memStatus == HIT_MS_NULL || tileRecv->memStatus == HIT_MS_NOMEM || hit_tileDims(*tileRecv)==0 ) return HIT_COM_NULL;
876 
877  if (selectRecvMode == HIT_COM_TILECOORDS) hit_tileSelect( &selectRecv, tileRecv, selectionRecv );
878  else if (selectRecvMode == HIT_COM_ARRAYCOORDS) hit_tileSelectArrayCoords( &selectRecv, tileRecv, selectionRecv);
879  else return HIT_COM_NULL;
880 
881  newCom.dataRecv = selectRecv.data;
882 
883  /* 6. COMMUNICATION TYPE */
884  newCom.commType = HIT_ALLTOALL;
885 
886  /* 7. TYPE and COUNT */
887  newCom.typeSend = baseType;
888  newCom.count = count;
889 
890  /* 8. RETURN */
891  return newCom;
892 }
893 
894 
895 
897  const void * tilePSend,
898  HitShape * selectionSend,
899  int selectSendMode,
900  const void * tilePRecv,
901  HitShape * selectionRecv,
902  int selectRecvMode,
903  HitType baseType){
904 
905 
906  /* 2. DECLARE NEW COMMUNICATION ISSUE */
907  HitCom newCom = HIT_COM_NULL;
908 
910  // @arturo Ago 2015: New allocP interface
911  // hit_malloc(newCom.alltoallv,sizeof(HitComAlltoallv),HitComAlltoallv*);
912  *(newCom.alltoallv) = HIT_COM_ALLTOALLV_NULL;
913 
914  /* 3. GET TOPOLOGY COMMUNICATOR */
915  //MPI_Comm comm = *((MPI_Comm *)lay.pTopology[0].lowLevel);
916  MPI_Comm comm = lay.pTopology[0]->comm;
917  newCom.comm = comm;
918 
919  /* 4. START FILLING UP ITS FIELDS */
920  // @arturo 2015/01/22
921  //newCom.myself = lay.topo.linearRank;
922  newCom.myself = hit_topSelfRankInternal( lay.topo );
923 
924 
925  /* 5. COMMUNICATION */
926  int numProcs = 1;
927  int i;
928  for(i=0;i<lay.topo.numDims;i++){
929  numProcs *= lay.numActives[i];
930  }
931 
932  hit_malloc(newCom.alltoallv->sendcnts, int, numProcs );
933  hit_malloc(newCom.alltoallv->sdispls, int, numProcs );
934  hit_malloc(newCom.alltoallv->recvcnts, int, numProcs );
935  hit_malloc(newCom.alltoallv->rdispls, int, numProcs );
936  // @arturo Ago 2015: New allocP interface
937  /*
938  hit_malloc(newCom.alltoallv->sendcnts,sizeof(int) * (size_t) numProcs,int*);
939  hit_malloc(newCom.alltoallv->sdispls ,sizeof(int) * (size_t) numProcs,int*);
940  hit_malloc(newCom.alltoallv->recvcnts,sizeof(int) * (size_t) numProcs,int*);
941  hit_malloc(newCom.alltoallv->rdispls ,sizeof(int) * (size_t) numProcs,int*);
942  */
943 
944 
945  const HitTile *tileSend = (const HitTile *)tilePSend;
946  newCom.dataSend = tileSend->data;
947 
948  const HitTile *tileRecv = (const HitTile *)tilePRecv;
949  newCom.dataRecv = tileRecv->data;
950 
951 
952 
953 
954  for(i=0; i<numProcs; i++){
955 
956  // Send
957  newCom.alltoallv->sendcnts[i] = hit_sigCard( hit_shapeSig(selectionSend[i],0));
958  // With hit_sigTileToArray we get the displacement of the first position,
959  // if HIT_COM_ARRAYCOORDS is selected the real displacement is calculated
960  // subtracting the displacement of the first element of the array.
961  if(selectSendMode == HIT_COM_TILECOORDS){
962  newCom.alltoallv->sdispls[i] = hit_sigTileToArray(hit_shapeSig(selectionSend[i],0),0);
963  } else {
964  newCom.alltoallv->sdispls[i] = hit_sigTileToArray(hit_shapeSig(selectionSend[i],0),0) -
965  hit_sigTileToArray(hit_shapeSig(tileSend->shape,0),0);
966  }
967 
968  // Recv
969  newCom.alltoallv->recvcnts[i] = hit_sigCard(hit_shapeSig(selectionRecv[i],0));
970  if(selectRecvMode == HIT_COM_TILECOORDS){
971  newCom.alltoallv->rdispls[i] = hit_sigTileToArray(hit_shapeSig(selectionRecv[i],0),0);
972  } else {
973  newCom.alltoallv->rdispls[i] = hit_sigTileToArray(hit_shapeSig(selectionRecv[i],0),0) -
974  hit_sigTileToArray(hit_shapeSig(tileRecv->shape,0),0);
975 
976  }
977  }
978 
979 
980 
981  /* 6. COMMUNICATION TYPE */
982  newCom.commType = HIT_ALLTOALLV;
983 
984  /* 7. TYPE */
985  newCom.typeSend = baseType;
986 
987  /* 8. RETURN */
988  return newCom;
989 
990 
991 }
992 
993 
994 
995 // @note @javfres Not used anymore in MatMult.
996 HitCom hit_comAllGathervInternal(HitLayout lay, const void * tilePSend, const void * tilePRecv, HitType baseType, const char * file, int line){
997 
998  HIT_NOT_USED(file);
999  HIT_NOT_USED(line);
1000 
1001  /* 2. DECLARE NEW COMMUNICATION ISSUE */
1002  HitCom newCom = HIT_COM_NULL;
1003 
1004  hit_malloc(newCom.alltoallv, HitComAlltoallv, 1);
1005  // @arturo Ago 2015: New allocP interface
1006  // hit_malloc(newCom.alltoallv,sizeof(HitComAlltoallv),HitComAlltoallv*);
1007  *(newCom.alltoallv) = HIT_COM_ALLTOALLV_NULL;
1008 
1009  /* 3. GET TOPOLOGY COMMUNICATOR */
1010  //MPI_Comm comm = *((MPI_Comm *)lay.pTopology[0].lowLevel);
1011  MPI_Comm comm = lay.pTopology[0]->comm;
1012  newCom.comm = comm;
1013 
1014  /* 4. START FILLING UP ITS FIELDS */
1015  // @arturo 2015/01/22
1016  //newCom.myself = lay.topo.linearRank;
1017  newCom.myself = hit_topSelfRankInternal( lay.topo );
1018 
1019 
1020  /* 5. COMMUNICATION */
1021  int numProcs = 1;
1022  int i;
1023  for(i=0;i<lay.topo.numDims;i++){
1024  numProcs *= lay.numActives[i];
1025  }
1026 
1027  hit_malloc(newCom.alltoallv->recvcnts, int, numProcs );
1028  hit_malloc(newCom.alltoallv->rdispls, int, numProcs );
1029  /*
1030  hit_malloc(newCom.alltoallv->recvcnts,sizeof(int) * (size_t) numProcs,int*);
1031  hit_malloc(newCom.alltoallv->rdispls ,sizeof(int) * (size_t) numProcs,int*);
1032  */
1033 
1034  const HitTile *tileSend = (const HitTile *)tilePSend;
1035  newCom.dataSend = tileSend->data;
1036 
1037  const HitTile *tileRecv = (const HitTile *)tilePRecv;
1038  newCom.dataRecv = tileRecv->data;
1039 
1040 // SEND
1041  HitShape myshape = hit_layShape(lay);
1042  newCom.count = hit_sigCard(hit_shapeSig(myshape,0));
1043 
1044 
1045 // RECV
1046  int acum_disp = 0;
1047  for(i=0; i<numProcs; i++){
1048 
1049  HitRanks ranks = HIT_RANKS_NULL;
1050  ranks.rank[0] = i;
1051 
1052  HitShape shape = hit_layShapeOther(lay, ranks);
1053 
1054  //printf("[%d] From %d: [%d-%d]\n",newCom.myself,i,hit_shapeSig(shape,0).begin,hit_shapeSig(shape,0).end);
1055 
1056  newCom.alltoallv->recvcnts[i] = hit_sigCard(hit_shapeSig(shape,0));
1057  newCom.alltoallv->rdispls[i] = acum_disp;
1058  acum_disp += hit_sigCard(hit_shapeSig(shape,0));
1059 
1060  //printf("[%d] From %d: Recv[counts:%d, disp:%d]\n",newCom.myself,i,newCom.alltoallv->recvcnts[i],newCom.alltoallv->rdispls[i]);
1061  }
1062 
1063 
1064  /* 6. COMMUNICATION TYPE */
1065  newCom.commType = HIT_ALLGATHERV;
1066 
1067  /* 7. TYPE */
1068  newCom.typeSend = baseType;
1069  newCom.typeRecv = baseType;
1070 
1071  /* 8. RETURN */
1072  return newCom;
1073 
1074 
1075 
1076 }
1077 
1078 
1079 
1080 
1081 
1082 
1083 
1084 HitCom hit_comSparseUpdate(HitLayout lay, const void * tileP, HitType baseType){
1085 
1086  const HitTile *tile = (const HitTile *)tileP;
1087 
1088  // Switch to select the appropriate function.
1089  // @arturo Aug 2015: Change in a name of internal function
1090  switch(hit_tileClass(*tile)){
1091  case HIT_GC_TILE:
1092  return hit_comSparseUpdateCSR(lay,tileP,baseType);
1093  case HIT_GB_TILE:
1094  return hit_comSparseUpdateBitmap(lay,tileP,baseType);
1095  default:
1096  hit_errInternal(__func__, "Unsupported shape type", "", __FILE__, __LINE__);
1097  return HIT_COM_NULL;
1098  }
1099 }
1100 
1101 
1102 //#define DEBUG_LAY_METIS
1103 
1104 #ifdef DEBUG_LAY_METIS
1105 #define debug(...) { if( hit_Rank == 0 ) { printf(__VA_ARGS__); fflush(stdout); }}
1106 #define debugall(...) { printf(__VA_ARGS__);}
1107 #endif
1108 
1109 HitCom hit_comSparseUpdateCSR(HitLayout lay, const void * tileP, HitType baseType){
1110 
1111 
1112  /* 1. DECLARE NEW COMMUNICATION ISSUE */
1113  HitCom newCom = HIT_COM_NULL;
1114  //@javfres 2015-10-05 Fix for inactive processes
1115  if(!lay.active) return newCom;
1116  // @arturo Ago 2015: New allocP interface
1117  // hit_malloc(newCom.alltoallv,sizeof(HitComAlltoallv),HitComAlltoallv*);
1118  hit_malloc(newCom.alltoallv, HitComAlltoallv, 1);
1119  *(newCom.alltoallv) = HIT_COM_ALLTOALLV_NULL;
1120  // @arturo Ago 2015: New allocP interface
1121  // hit_malloc(newCom.alltoallv->sparse,sizeof(HitComSparse),HitComSparse*);
1122  hit_malloc(newCom.alltoallv->sparse, HitComSparse, 1);
1123  *(newCom.alltoallv->sparse) = HIT_COM_SPARSE_NULL;
1124 
1125 
1126  /* 2. GET TOPOLOGY COMMUNICATOR */
1127  //MPI_Comm comm = *((MPI_Comm *)lay.pTopology[0].lowLevel);
1128  MPI_Comm comm = lay.pTopology[0]->comm;
1129  newCom.comm = comm;
1130 
1131  /* 3. START FILLING UP ITS FIELDS */
1132  // @arturo 2015/01/22
1133  //newCom.myself = lay.topo.linearRank;
1134  newCom.myself = hit_topSelfRankInternal( lay.topo );
1135 
1136  /* 4. COMMUNICATION */
1137  int numProcs = 1;
1138  int i;
1139  for(i=0;i<lay.topo.numDims;i++){
1140  numProcs *= lay.numActives[i];
1141  }
1142 
1143 
1144  /* 6. HIT TILE */
1145  const HitTile *tile = (const HitTile *)tileP;
1146 
1147 
1148  HitShape ext_shape = tile->shape;
1149 
1150 
1151  /* 5. INIT THE SEND AND RECV ARRAYS */
1152 
1153  // Shape is the original shape.
1154  HitShape shape = lay.origShape;
1155 
1156  // Part is the Metis output.
1157  int * part = lay.info.layoutList.assignedGroups;
1158 
1159  // Local is an array with the local + neighbor vertices.
1160  int * local;
1161 
1162  // @arturo Ago 2015: New allocP interface
1163  // hit_malloc(local,sizeof(int) * (size_t) hit_cShapeNvertices(shape), int*);
1164  hit_malloc(local, int, hit_cShapeNvertices(shape) );
1165 
1166  for(i=0;i<hit_cShapeNvertices(shape);i++){
1167  if(part[i] == lay.group){
1168  local[i] = 1;
1169  } else {
1170  local[i] = 0;
1171  }
1172  }
1173 
1174  // Set as local the neighbor vertices
1175  for(i=hit_cShapeNvertices(lay.shape);i<hit_cShapeNvertices(ext_shape);i++){
1176 
1177  int gvertex = hit_cShapeVertexToGlobal(ext_shape,i);
1178  int overtex = hit_cShapeVertexToLocal(shape,gvertex);
1179  local[overtex] = 1;
1180  }
1181 
1182 
1183  // 5.1 Calculate the recv
1184  int * recvcnts;
1185  int * recvdispls;
1186  int * auxcounts;
1187 
1188  hit_calloc(recvcnts, int, numProcs);
1189  hit_calloc(recvdispls,int, numProcs);
1190  hit_calloc(auxcounts, int, numProcs);
1191  // @arturo Ago 2015: New allocP interface
1192  /*
1193  hit_calloc(recvcnts, sizeof(int), (size_t) numProcs,int*);
1194  hit_calloc(recvdispls,sizeof(int), (size_t) numProcs,int*);
1195  hit_calloc(auxcounts, sizeof(int), (size_t) numProcs,int*);
1196  */
1197 
1198  int nrecv = 0;
1199  for(i=0; i<hit_cShapeNvertices(shape); i++){
1200 
1201  if( local[i] && part[i] != lay.group){
1202  recvcnts[part[i]] ++;
1203  nrecv ++;
1204  }
1205  }
1206 
1207  recvdispls[0] = 0;
1208  for(i=1;i<numProcs;i++){
1209  recvdispls[i] = recvdispls[i-1] + recvcnts[i-1];
1210  }
1211 
1212 
1213  int * recv;
1214  hit_malloc(recv, int, nrecv);
1215  // @arturo Ago 2015: New allocP interface
1216  // hit_malloc(recv,sizeof(int) * (size_t) nrecv,int*);
1217 
1218  for(i=0; i<hit_cShapeNvertices(shape); i++){
1219 
1220  if( local[i] && part[i] != lay.group){
1221  int p = part[i];
1222  recv[ recvdispls[p] + auxcounts[p] ] = hit_cShapeNameList(shape,0).names[i];
1223  auxcounts[p]++;
1224  }
1225  }
1226 
1227 #ifdef DEBUG_LAY_METIS
1228  { sleep(hit_Rank); debugall("I recv: [");
1229  int p, i;
1230  for(p=0;p<numProcs;p++){
1231  debugall("p%d(%d): ",p,recvcnts[p]);
1232  for(i=0;i<recvcnts[p];i++){
1233  debugall("%d ",recv[ recvdispls[p] + i ]);
1234  }
1235  } debugall("]\n");}
1236 #endif
1237 
1238 
1239  // 6.2 Calculate the send
1240  int * sendcnts;
1241  int * senddispls;
1242 
1243  hit_malloc(sendcnts, int, numProcs);
1244  hit_malloc(senddispls, int, numProcs);
1245  // @arturo Ago 2015: New allocP interface
1246  /*
1247  hit_malloc(sendcnts,sizeof(int) * (size_t) numProcs,int*);
1248  hit_malloc(senddispls,sizeof(int) * (size_t) numProcs,int*);
1249  */
1250 
1251  int ok = MPI_Alltoall(recvcnts, 1, MPI_INT, sendcnts, 1, MPI_INT, comm);
1252  hit_mpiTestError(ok,"Failed while sending the recvcnts");
1253 
1254  senddispls[0] = 0;
1255  for(i=1;i<numProcs;i++){
1256  senddispls[i] = senddispls[i-1] + sendcnts[i-1];
1257  }
1258 
1259 
1260  int * send;
1261  hit_malloc(send, int, senddispls[numProcs-1] + sendcnts[numProcs-1] );
1262  // @arturo Ago 2015: New allocP interface
1263  // hit_malloc(send,sizeof(int) * (size_t) ( senddispls[numProcs-1] + sendcnts[numProcs-1] ),int*);
1264 
1265 
1266  MPI_Alltoallv(recv, recvcnts, recvdispls,
1267  MPI_INT, send, sendcnts,
1268  senddispls, MPI_INT, comm);
1269  hit_mpiTestError(ok,"Failed while sending the recv list of vertices.");
1270 
1271 
1272 #ifdef DEBUG_LAY_METIS
1273  { sleep(hit_Rank); debugall("I send: [");
1274  int p, i;
1275  for(p=0;p<numProcs;p++){
1276  debugall("p%d(%d): ",p,sendcnts[p]);
1277  for(i=0;i<sendcnts[p];i++){
1278  debugall("%d ",send[ senddispls[p] + i ]);
1279  }
1280  } debugall("]\n");}
1281 #endif
1282 
1283 
1284  newCom.alltoallv->sendcnts = sendcnts;
1285  newCom.alltoallv->sdispls = senddispls;
1286  newCom.alltoallv->recvcnts = recvcnts;
1287  newCom.alltoallv->rdispls = recvdispls;
1288  newCom.alltoallv->sparse->send = send;
1289  newCom.alltoallv->sparse->recv = recv;
1290  newCom.alltoallv->sparse->nsend = newCom.alltoallv->sdispls[numProcs-1] + newCom.alltoallv->sendcnts[numProcs-1];
1291  newCom.alltoallv->sparse->nrecv = newCom.alltoallv->rdispls[numProcs-1] + newCom.alltoallv->recvcnts[numProcs-1];
1292 
1293 
1294  hit_vmalloc(newCom.dataSend, tile->baseExtent * (size_t) newCom.alltoallv->sparse->nsend );
1295  hit_vmalloc(newCom.dataRecv, tile->baseExtent * (size_t) newCom.alltoallv->sparse->nrecv );
1296  // @arturo Ago 2015: New allocP interface
1297  /*
1298  hit_malloc(newCom.dataSend,tile->baseExtent * (size_t) newCom.alltoallv->sparse->nsend,void*);
1299  hit_malloc(newCom.dataRecv,tile->baseExtent * (size_t) newCom.alltoallv->sparse->nrecv,void*);
1300  */
1301 
1302  /* 7. FILL THE SPARSE RELATED FIELDS */
1303  newCom.alltoallv->sparse->shape = tile->shape;
1304  newCom.alltoallv->sparse->originData = tile->dataVertices;
1305 
1306 
1307  /* 8. COMMUNICATION TYPE */
1308  newCom.commType = HIT_SPARSE_UPDATE;
1309 
1310  /* 9. TYPE */
1311  newCom.typeSend = baseType;
1312  newCom.alltoallv->sparse->baseExtent = tile->baseExtent;
1313 
1314  // 10. Release resources.
1315  free(auxcounts);
1316  free(local);
1317 
1318 
1319  /* 11. RETURN */
1320  return newCom;
1321 }
1322 
1323 
1324 HitCom hit_comSparseUpdateBitmap(HitLayout lay, const void * tileP, HitType baseType){
1325 
1326  /* 1. DECLARE NEW COMMUNICATION ISSUE */
1327  HitCom newCom = HIT_COM_NULL;
1328  //@javfres 2015-10-05 Fix for inactive processes
1329  if(!lay.active) return newCom;
1330  // @arturo Ago 2015: New allocP interface
1331  // hit_malloc(newCom.alltoallv,sizeof(HitComAlltoallv),HitComAlltoallv*);
1332  hit_malloc(newCom.alltoallv, HitComAlltoallv, 1);
1333  *(newCom.alltoallv) = HIT_COM_ALLTOALLV_NULL;
1334  // @arturo Ago 2015: New allocP interface
1335  // hit_malloc(newCom.alltoallv->sparse,sizeof(HitComSparse),HitComSparse*);
1336  hit_malloc(newCom.alltoallv->sparse, HitComSparse, 1);
1337  *(newCom.alltoallv->sparse) = HIT_COM_SPARSE_NULL;
1338 
1339 
1340  /* 2. GET TOPOLOGY COMMUNICATOR */
1341  //MPI_Comm comm = *((MPI_Comm *)lay.pTopology[0].lowLevel);
1342  MPI_Comm comm = lay.pTopology[0]->comm;
1343  newCom.comm = comm;
1344 
1345  /* 3. START FILLING UP ITS FIELDS */
1346  // @arturo 2015/01/22
1347  //newCom.myself = lay.topo.linearRank;
1348  newCom.myself = hit_topSelfRankInternal( lay.topo );
1349 
1350  /* 4. COMMUNICATION */
1351  int numProcs = 1;
1352  int i;
1353  for(i=0;i<lay.topo.numDims;i++){
1354  numProcs *= lay.numActives[i];
1355  }
1356 
1357 
1358  /* 6. HIT TILE */
1359  const HitTile *tile = (const HitTile *)tileP;
1360 
1361 
1362  HitShape ext_shape = tile->shape;
1363 
1364 
1365  /* 5. INIT THE SEND AND RECV ARRAYS */
1366 
1367  // Shape is the original shape.
1368  HitShape shape = lay.origShape;
1369 
1370  // Part is the Metis output.
1371  int * part = lay.info.layoutList.assignedGroups;
1372 
1373  // Local is an array with the local + neighbor vertices.
1374  int * local;
1375 
1376  // @arturo Ago 2015: New allocP interface
1377  // hit_malloc(local,sizeof(int) * (size_t) hit_bShapeNvertices(shape),int*);
1378  hit_malloc(local, int, hit_bShapeNvertices(shape));
1379 
1380  for(i=0;i<hit_bShapeNvertices(shape);i++){
1381  if(part[i] == lay.group){
1382  local[i] = 1;
1383  } else {
1384  local[i] = 0;
1385  }
1386  }
1387 
1388  // Set as local the neighbor vertices
1389  for(i=hit_bShapeNvertices(lay.shape);i<hit_bShapeNvertices(ext_shape);i++){
1390 
1391  int gvertex = hit_bShapeVertexToGlobal(ext_shape,i);
1392  int overtex = hit_bShapeVertexToLocal(shape,gvertex);
1393  local[overtex] = 1;
1394  }
1395 
1396 
1397  // 5.1 Calculate the recv
1398  int * recvcnts;
1399  int * recvdispls;
1400  int * auxcounts;
1401 
1402  hit_calloc(recvcnts, int, numProcs);
1403  hit_calloc(recvdispls,int, numProcs);
1404  hit_calloc(auxcounts, int, numProcs);
1405  // @arturo Ago 2015: New allocP interface
1406  /*
1407  hit_calloc(recvcnts, sizeof(int), (size_t) numProcs,int*);
1408  hit_calloc(recvdispls,sizeof(int), (size_t) numProcs,int*);
1409  hit_calloc(auxcounts, sizeof(int), (size_t) numProcs,int*);
1410  */
1411 
1412 
1413  int nrecv = 0;
1414  for(i=0; i<hit_bShapeNvertices(shape); i++){
1415 
1416  if( local[i] && part[i] != lay.group){
1417  recvcnts[part[i]] ++;
1418  nrecv ++;
1419  }
1420  }
1421 
1422  recvdispls[0] = 0;
1423  for(i=1;i<numProcs;i++){
1424  recvdispls[i] = recvdispls[i-1] + recvcnts[i-1];
1425  }
1426 
1427 
1428  int * recv;
1429  // @arturo Ago 2015: New allocP interface
1430  // hit_malloc(recv,sizeof(int) * (size_t) nrecv, int*);
1431  hit_malloc(recv, int, nrecv);
1432 
1433  for(i=0; i<hit_bShapeNvertices(shape); i++){
1434 
1435  if( local[i] && part[i] != lay.group){
1436  int p = part[i];
1437  recv[ recvdispls[p] + auxcounts[p] ] = hit_bShapeNameList(shape,0).names[i];
1438  auxcounts[p]++;
1439  }
1440  }
1441 
1442 #ifdef DEBUG_LAY_METIS
1443  { sleep(hit_Rank); debugall("I recv: [");
1444  int p, i;
1445  for(p=0;p<numProcs;p++){
1446  debugall("p%d(%d): ",p,recvcnts[p]);
1447  for(i=0;i<recvcnts[p];i++){
1448  debugall("%d ",recv[ recvdispls[p] + i ]);
1449  }
1450  } debugall("]\n");}
1451 #endif
1452 
1453 
1454  // 6.2 Calculate the send
1455  int * sendcnts;
1456  int * senddispls;
1457 
1458  hit_malloc(sendcnts, int, numProcs);
1459  hit_malloc(senddispls, int, numProcs);
1460  /*
1461  hit_malloc(sendcnts,sizeof(int) * (size_t) numProcs,int*);
1462  hit_malloc(senddispls,sizeof(int) * (size_t) numProcs,int*);
1463  */
1464 
1465  int ok = MPI_Alltoall(recvcnts, 1, MPI_INT, sendcnts, 1, MPI_INT, comm);
1466  hit_mpiTestError(ok,"Failed while sending the recvcnts");
1467 
1468  senddispls[0] = 0;
1469  for(i=1;i<numProcs;i++){
1470  senddispls[i] = senddispls[i-1] + sendcnts[i-1];
1471  }
1472 
1473 
1474  int * send;
1475  // @arturo Ago 2015: New allocP interface
1476  // hit_malloc(send,sizeof(int) * (size_t) ( senddispls[numProcs-1] + sendcnts[numProcs-1] ),int*);
1477  hit_malloc(send, int, senddispls[numProcs-1] + sendcnts[numProcs-1] );
1478 
1479 
1480  MPI_Alltoallv(recv, recvcnts, recvdispls,
1481  MPI_INT, send, sendcnts,
1482  senddispls, MPI_INT, comm);
1483  hit_mpiTestError(ok,"Failed while sending the recv list of vertices.");
1484 
1485 
1486 
1487 #ifdef DEBUG_LAY_METIS
1488  { sleep(hit_Rank); debugall("I send: [");
1489  int p, i;
1490  for(p=0;p<numProcs;p++){
1491  debugall("p%d(%d): ",p,sendcnts[p]);
1492  for(i=0;i<sendcnts[p];i++){
1493  debugall("%d ",send[ senddispls[p] + i ]);
1494  }
1495  } debugall("]\n");}
1496 #endif
1497 
1498 
1499 
1500  newCom.alltoallv->sendcnts = sendcnts;
1501  newCom.alltoallv->sdispls = senddispls;
1502  newCom.alltoallv->recvcnts = recvcnts;
1503  newCom.alltoallv->rdispls = recvdispls;
1504  newCom.alltoallv->sparse->send = send;
1505  newCom.alltoallv->sparse->recv = recv;
1506  newCom.alltoallv->sparse->nsend = newCom.alltoallv->sdispls[numProcs-1] + newCom.alltoallv->sendcnts[numProcs-1];
1507  newCom.alltoallv->sparse->nrecv = newCom.alltoallv->rdispls[numProcs-1] + newCom.alltoallv->recvcnts[numProcs-1];
1508 
1509 
1510  hit_vmalloc(newCom.dataSend, tile->baseExtent * (size_t) newCom.alltoallv->sparse->nsend);
1511  hit_vmalloc(newCom.dataRecv, tile->baseExtent * (size_t) newCom.alltoallv->sparse->nrecv);
1512  // @arturo Ago 2015: New allocP interface
1513  /*
1514  hit_malloc(newCom.dataSend,tile->baseExtent * (size_t) newCom.alltoallv->sparse->nsend,void*);
1515  hit_malloc(newCom.dataRecv,tile->baseExtent * (size_t) newCom.alltoallv->sparse->nrecv,void*);
1516  */
1517 
1518  /* 7. FILL THE SPARSE RELATED FIELDS */
1519  newCom.alltoallv->sparse->shape = tile->shape;
1520  newCom.alltoallv->sparse->originData = tile->dataVertices;
1521 
1522  /* 8. COMMUNICATION TYPE */
1523  newCom.commType = HIT_SPARSE_UPDATE;
1524 
1525  /* 9. TYPE */
1526  newCom.typeSend = baseType;
1527  newCom.alltoallv->sparse->baseExtent = tile->baseExtent;
1528 
1529  // 10. Release resources.
1530  free(auxcounts);
1531  free(local);
1532 
1533 
1534  /* 11. RETURN */
1535  return newCom;
1536 }
1537 
1538 
1539 
1540 
1541 
1542 
1543 HitCom hit_comSparseScatterInternal(HitLayout lay, const void * tilePSend, const void * tilePRecv, HitType baseType, const char * file, int line){
1544 
1545  if(!lay.active) return HIT_COM_NULL;
1546 
1547 
1548  /* 1. HIT TILE */
1549  const HitTile *tileSend = (const HitTile *) tilePSend;
1550  const HitTile *tileRecv = (const HitTile *) tilePRecv;
1551 
1552 
1553  // @arturo 2015/01/22
1554  //if(lay.topo.linearRank == 0 && hit_tileType(*tileSend) != hit_tileType(*tileRecv)){
1555  // @arturo Aug 2015: Change in a name of internal function
1556  //if( hit_topSelfRankInternal( lay.topo ) == 0 && hit_tileType(*tileSend) != hit_tileType(*tileRecv)){
1557  if( hit_topSelfRankInternal( lay.topo ) == 0 && hit_tileClass(*tileSend) != hit_tileClass(*tileRecv)){
1558  hit_errInternal(__func__, "Different HitTile classes", "", file, line);
1559  }
1560 
1561  //if(hit_tileType(*tileRecv) != HIT_GC_TILE && hit_tileType(*tileRecv) != HIT_GB_TILE){
1562  if(hit_tileClass(*tileRecv) != HIT_GC_TILE && hit_tileClass(*tileRecv) != HIT_GB_TILE){
1563  hit_errInternal(__func__, "Unsupported HitTile class", "", __FILE__, __LINE__);
1564  }
1565 
1566  /* 2. DECLARE NEW COMMUNICATION ISSUE */
1567  HitCom newCom = HIT_COM_NULL;
1568  // @arturo Ago 2015: New allocP interface
1569  // hit_malloc(newCom.alltoallv,sizeof(HitComAlltoallv),HitComAlltoallv*);
1570  hit_malloc(newCom.alltoallv, HitComAlltoallv, 1);
1571  *(newCom.alltoallv) = HIT_COM_ALLTOALLV_NULL;
1572  // @arturo Ago 2015: New allocP interface
1573  // hit_malloc(newCom.alltoallv->sparse,sizeof(HitComSparse),HitComSparse*);
1574  hit_malloc(newCom.alltoallv->sparse, HitComSparse, 1);
1575  *(newCom.alltoallv->sparse) = HIT_COM_SPARSE_NULL;
1576 
1577 
1578  /* 3. GET TOPOLOGY COMMUNICATOR */
1579  //MPI_Comm comm = *((MPI_Comm *)lay.pTopology[0].lowLevel);
1580  MPI_Comm comm = lay.pTopology[0]->comm;
1581  newCom.comm = comm;
1582 
1583  /* 4. START FILLING UP ITS FIELDS */
1584  // @arturo 2015/01/22
1585  //newCom.myself = lay.topo.linearRank;
1586  newCom.myself = hit_topSelfRankInternal( lay.topo );
1587 
1588  /* 5. COMMUNICATION */
1589  int numProcs = 1;
1590  int i;
1591  for(i=0;i<lay.topo.numDims;i++){
1592  numProcs *= lay.numActives[i];
1593  }
1594 
1595  /* 5. INIT THE SEND AND RECV ARRAYS */
1596  // @arturo Aug 2015: Change in a name of internal function
1597  //if(hit_tileType(*tileRecv) == HIT_GC_TILE){
1598  if(hit_tileClass(*tileRecv) == HIT_GC_TILE){
1599  newCom.count = hit_cShapeNvertices(lay.shape);
1600  } else {
1601  newCom.count = hit_bShapeNvertices(lay.shape);
1602  }
1603 
1604 
1605  if(newCom.myself == 0){
1606 
1607  hit_calloc(newCom.alltoallv->sendcnts, int, numProcs);
1608  hit_malloc(newCom.alltoallv->sdispls, int, numProcs);
1609  // @arturo Ago 2015: New allocP interface
1610  /*
1611  hit_calloc(newCom.alltoallv->sendcnts,(size_t) numProcs, sizeof(int),int*);
1612  hit_malloc(newCom.alltoallv->sdispls,sizeof(int) * (size_t) numProcs,int*);
1613  */
1614 
1615  for(i=0;i<lay.info.layoutList.numElementsTotal;i++){
1616  newCom.alltoallv->sendcnts[ lay.info.layoutList.assignedGroups[i] ]++;
1617  }
1618 
1619  newCom.alltoallv->sdispls[0] = 0;
1620  for(i=1;i<numProcs;i++){
1621  newCom.alltoallv->sdispls[i] = newCom.alltoallv->sdispls[i-1] + newCom.alltoallv->sendcnts[i-1];
1622  }
1623 
1624  int total_size = newCom.alltoallv->sdispls[numProcs-1] + newCom.alltoallv->sendcnts[numProcs-1];
1625  // @arturo Ago 2015: New allocP interface
1626  // hit_malloc(newCom.dataSend,tileSend->baseExtent * (size_t) total_size,void*);
1627  hit_vmalloc(newCom.dataSend, tileSend->baseExtent * (size_t) total_size);
1628 
1629  newCom.alltoallv->sparse->originData = tileSend->dataVertices;
1632  }
1633  newCom.dataRecv = tileRecv->dataVertices;
1634 
1635 
1636  /* 8. COMMUNICATION TYPE */
1637  newCom.commType = HIT_SPARSE_SCATTER;
1638 
1639  /* 9. TYPE */
1640  newCom.typeSend = baseType;
1641  newCom.typeRecv = baseType;
1642 
1643  newCom.alltoallv->sparse->baseExtent = tileSend->baseExtent;
1644 
1645  /* 10. RETURN */
1646  return newCom;
1647 
1648 }
1649 
1650 
1652  const void * tilePSend,
1653  int * count,
1654  const void * tilePRecv,
1655  HitType baseType){
1656 
1657 
1658  /* 1. CHECK THE LAYOUT TYPE */
1659  if(lay.type != HIT_LAYOUT_CONTIGUOUS){
1660  hit_errInternal("hit_comAlldistribute", "This communication works only with a contiguous layout list", "", __FILE__, __LINE__);
1661  }
1662 
1663  /* 2. DECLARE NEW COMMUNICATION ISSUE */
1664  HitCom newCom = HIT_COM_NULL;
1665  // @arturo Ago 2015: New allocP interface
1666  // hit_malloc(newCom.alltoallv,sizeof(HitComAlltoallv),HitComAlltoallv*);
1667  hit_malloc(newCom.alltoallv, HitComAlltoallv, 1);
1668  *(newCom.alltoallv) = HIT_COM_ALLTOALLV_NULL;
1669 
1670  /* 3. GET TOPOLOGY COMMUNICATOR */
1671  //MPI_Comm comm = *((MPI_Comm *)lay.pTopology[0].lowLevel);
1672  MPI_Comm comm = lay.pTopology[0]->comm;
1673  newCom.comm = comm;
1674 
1675  /* 4. START FILLING UP ITS FIELDS */
1676  // @arturo 2015/01/22
1677  //newCom.myself = lay.topo.linearRank;
1678  newCom.myself = hit_topSelfRankInternal( lay.topo );
1679  newCom.dataSend = ((const HitTile *)tilePSend)->data;
1680  newCom.dataRecv = ((const HitTile *)tilePRecv)->data;
1681 
1682  /* 5. CALCULATE NUMBER OF SETS */
1683  /*
1684  int i;
1685  int nSets = 1;
1686  for(i=0;i<hit_shapeDims(lay.origShape);i++){
1687  nSets *= hit_sigCard(hit_shapeSig(lay.origShape,i));
1688  }
1689  */
1690  int nSets = hit_shapeCard( lay.origShape );
1691 
1692  /* 6. COMMUNICATION OF PART SIZES */
1693  int i;
1694  int numProcs = 1;
1695  for(i=0;i<lay.topo.numDims;i++){
1696  numProcs *= lay.numActives[i];
1697  }
1698 
1699  hit_malloc(newCom.alltoallv->sendcnts, int, numProcs);
1700  hit_malloc(newCom.alltoallv->sdispls , int, numProcs);
1701  hit_malloc(newCom.alltoallv->recvcnts, int, numProcs);
1702  hit_malloc(newCom.alltoallv->rdispls , int, numProcs);
1703  // @arturo Ago 2015: New allocP interface
1704  /*
1705  hit_malloc(newCom.alltoallv->sendcnts,sizeof(int) * (size_t) numProcs,int*);
1706  hit_malloc(newCom.alltoallv->sdispls ,sizeof(int) * (size_t) numProcs,int*);
1707  hit_malloc(newCom.alltoallv->recvcnts,sizeof(int) * (size_t) numProcs,int*);
1708  hit_malloc(newCom.alltoallv->rdispls ,sizeof(int) * (size_t) numProcs,int*);
1709  */
1710 
1711  newCom.alltoallv->sdispls[0] = 0;
1712  newCom.alltoallv->sendcnts[0] = 0;
1713  int g,g_ant;
1714  g = g_ant = 0;
1715  for(i=0;i<nSets;i++){
1716 
1717  g = hit_lgr_elementGroup(lay,i);
1718 
1719  if(g != g_ant){
1720  newCom.alltoallv->sdispls[g] = newCom.alltoallv->sdispls[g-1] + newCom.alltoallv->sendcnts[g-1];
1721  newCom.alltoallv->sendcnts[g] = 0;
1722  g_ant = g;
1723  }
1724  newCom.alltoallv->sendcnts[g] += count[i];
1725 
1726  }
1727 
1728  // @arturo: 2012-Jun, BUG CORRECTED (memory overrun detected with valgrind)
1729  //while( g < numProcs ){
1730  while( g < numProcs-1 ){
1731  g++;
1732  newCom.alltoallv->sendcnts[g] = 0;
1733  }
1734 
1735  int ok = MPI_Alltoall(newCom.alltoallv->sendcnts, 1, MPI_INT,newCom.alltoallv->recvcnts, 1, MPI_INT,comm);
1736  hit_mpiTestError(ok,"hit_comAllDistribute: Failed Alltoall");
1737 
1738  /* 7. COMPUTE DISPLACEMENTS AND SIZE OF RECEIVING BUFFER */
1739  newCom.alltoallv->rdispls[0] = 0;
1740  newCom.count = newCom.alltoallv->recvcnts[0];
1741  for(i=1; i<numProcs; i++) {
1742  newCom.alltoallv->rdispls[i] = newCom.alltoallv->rdispls[i-1] + newCom.alltoallv->recvcnts[i-1];
1743  newCom.count += newCom.alltoallv->recvcnts[i];
1744  }
1745 
1746  /* DEBUG ASSERT: CHECK THE SIZE OF THE RECEIVING BUFFER */
1747  if( newCom.count < hit_tileCard( *(const HitTile *)tilePRecv ) ){
1748  hit_errInternal("hit_comAlldistribute", "The receiving buffer does not have enough cardinality", "", __FILE__, __LINE__);
1749  }
1750 
1751  /* 8. COMMUNICATION TYPE */
1752  newCom.commType = HIT_ALLTOALLV; //HIT_ALLDISTRIBUTE;
1753 
1754  /* 8. TYPE */
1755  newCom.typeSend = baseType;
1756 
1757  /* 9. RETURN */
1758  return newCom;
1759 
1760 }
1761 
1762 
1763 
1764 HitCom hit_comSparseScatterRowsInternal(HitLayout lay, const void * tilePSend, const void * tilePRecv, HitType baseType, const char * file, int line){
1765 
1766  /* 1. HIT TILE */
1767  const HitTile *tileSend = (const HitTile *) tilePSend;
1768  const HitTile *tileRecv = (const HitTile *) tilePRecv;
1769 
1770 
1771  // @arturo 2015/01/22
1772  //if(lay.topo.linearRank == 0 && hit_tileType(*tileSend) != hit_tileType(*tileRecv)){
1773  // @arturo Aug 2015: Change in a name of internal function
1774  //if( hit_topSelfRankInternal( lay.topo ) == 0 && hit_tileType(*tileSend) != hit_tileType(*tileRecv)){
1775  if( hit_topSelfRankInternal( lay.topo ) == 0 && hit_tileClass(*tileSend) != hit_tileClass(*tileRecv)){
1776  hit_errInternal(__func__, "Different HitTile classes", "", file, line);
1777  }
1778 
1779  // @arturo Aug 2015: Change in a name of internal function
1780  //if(hit_tileType(*tileRecv) != HIT_MC_TILE ){
1781  if(hit_tileClass(*tileRecv) != HIT_MC_TILE ){
1782  hit_errInternal(__func__, "Unsupported HitTile class", "", __FILE__, __LINE__);
1783  }
1784 
1785  /* 2. DECLARE NEW COMMUNICATION ISSUE */
1786  HitCom newCom = HIT_COM_NULL;
1787  // @arturo Ago 2015: New allocP interface
1788  // hit_malloc(newCom.alltoallv,sizeof(HitComAlltoallv),HitComAlltoallv*);
1789  hit_malloc(newCom.alltoallv, HitComAlltoallv, 1);
1790  *(newCom.alltoallv) = HIT_COM_ALLTOALLV_NULL;
1791  // @arturo Ago 2015: New allocP interface
1792  // hit_malloc(newCom.alltoallv->sparse,sizeof(HitComSparse),HitComSparse*);
1793  hit_malloc(newCom.alltoallv->sparse, HitComSparse, 1);
1794  *(newCom.alltoallv->sparse) = HIT_COM_SPARSE_NULL;
1795 
1796 
1797  /* 3. GET TOPOLOGY COMMUNICATOR */
1798  //MPI_Comm comm = *((MPI_Comm *)lay.pTopology[0].lowLevel);
1799  MPI_Comm comm = lay.pTopology[0]->comm;
1800  newCom.comm = comm;
1801 
1802  /* 4. START FILLING UP ITS FIELDS */
1803  // @arturo 2015/01/22
1804  //newCom.myself = lay.topo.linearRank;
1805  newCom.myself = hit_topSelfRankInternal( lay.topo );
1806 
1807  /* 5. COMMUNICATION */
1808  int numProcs = 1;
1809  int i;
1810  for(i=0;i<lay.topo.numDims;i++){
1811  numProcs *= lay.numActives[i];
1812  }
1813 
1814  /* 5. INIT THE SEND AND RECV ARRAYS */
1815  newCom.count = hit_cShapeNZElems(lay.shape);
1816 
1817 
1818  if(newCom.myself == 0){
1819 
1820  hit_calloc(newCom.alltoallv->sendcnts, int, numProcs);
1821  hit_malloc(newCom.alltoallv->sdispls, int, numProcs);
1822  /*
1823  hit_calloc(newCom.alltoallv->sendcnts,(size_t) numProcs, sizeof(int),int*);
1824  hit_malloc(newCom.alltoallv->sdispls,sizeof(int) * (size_t) numProcs,int*);
1825  */
1826 
1827  for(i=0;i<lay.info.layoutList.numElementsTotal;i++){
1828  int p = lay.info.layoutList.assignedGroups[i];
1829 
1830  newCom.alltoallv->sendcnts[ p ] += hit_cShapeNColsRow(lay.origShape,i);
1831 
1832  //printf("Send row[%d] with %d cols to %d\n",i,hit_cShapeNColsRow(lay.origShape,i),p);
1833  }
1834 
1835  newCom.alltoallv->sdispls[0] = 0;
1836  for(i=1;i<numProcs;i++){
1837  newCom.alltoallv->sdispls[i] = newCom.alltoallv->sdispls[i-1] + newCom.alltoallv->sendcnts[i-1];
1838  }
1839 
1840  int total_size = newCom.alltoallv->sdispls[numProcs-1] + newCom.alltoallv->sendcnts[numProcs-1];
1841  // @arturo Ago 2015: New allocP interface
1842  // hit_malloc(newCom.dataSend,tileSend->baseExtent * (size_t) total_size,void*);
1843  hit_vmalloc(newCom.dataSend, tileSend->baseExtent * (size_t) total_size );
1844 
1845  newCom.alltoallv->sparse->originData = tileSend->data;
1848 
1849  newCom.alltoallv->sparse->rows = hit_cShapeXadj(lay.origShape);
1850  }
1851  newCom.dataRecv = tileRecv->data;
1852 
1853 
1854  /* 8. COMMUNICATION TYPE */
1856 
1857  /* 9. TYPE */
1858  newCom.typeSend = baseType;
1859  newCom.typeRecv = baseType;
1860 
1861  newCom.alltoallv->sparse->baseExtent = tileSend->baseExtent;
1862 
1863  /* 10. RETURN */
1864  return newCom;
1865 
1866 }
1867 
1868 
1869 #ifdef XXXXXXX
1871  const void * tilePSend, int * count,
1872  const void * tilePRecv, HitType baseType){
1873 
1874 
1875  if(lay.type != HIT_LAYOUT_CONTIGUOUS){
1876  hit_errInternal("hit_comAlldistribute", "This communication works only with a contiguous layout list", "", "", -1);
1877  }
1878 
1879  /* 2. DECLARE NEW COMMUNICATION ISSUE */
1880  HitCom newCom = HIT_COM_NULL;
1881  // @arturo Ago 2015: New allocP interface
1882  // hit_malloc(newCom.alltoallv,sizeof(HitComAlltoallv));
1883  hit_malloc(newCom.alltoallv, HitComAlltoallv, 1);
1884  *(newCom.alltoallv) = HIT_COM_ALLTOALLV_NULL;
1885 
1886 
1887  /* 3. GET TOPOLOGY COMMUNICATOR */
1888  MPI_Comm comm = *((MPI_Comm *)lay.pTopology[0].lowLevel);
1889  newCom.comm = comm;
1890 
1891  /* 4. START FILLING UP ITS FIELDS */
1892  // @arturo 2015/01/22
1893  //newCom.myself = lay.topo.linearRank;
1894  newCom.myself = hit_topSelfRankInternal( lay.topo );
1895 
1896  int i;
1897 
1898  /* CALCULATE NUMBER OF SETS */
1899  int nSets = 1;
1900  for(i=0;i<hit_shapeDims(lay.origShape);i++){
1901  nSets *= hit_sigCard(hit_shapeSig(lay.origShape,i));
1902  }
1903 
1904  /* 5. COMMUNICATION */
1905  int numProcs = 1;
1906  for(i=0;i<lay.topo.numDims;i++){
1907  numProcs *= lay.numActives[i];
1908  }
1909 
1910 
1911  hit_malloc(newCom.alltoallv->sendcnts, int, numProcs);
1912  hit_malloc(newCom.alltoallv->sdispls , int, numProcs);
1913  hit_malloc(newCom.alltoallv->recvcnts, int, numProcs);
1914  hit_malloc(newCom.alltoallv->rdispls , int, numProcs);
1915  // @arturo Ago 2015: New allocP interface
1916  /*
1917  hit_malloc(newCom.alltoallv->sendcnts,sizeof(int) * (size_t) numProcs);
1918  hit_malloc(newCom.alltoallv->sdispls ,sizeof(int) * (size_t) numProcs);
1919  hit_malloc(newCom.alltoallv->recvcnts,sizeof(int) * (size_t) numProcs);
1920  hit_malloc(newCom.alltoallv->rdispls ,sizeof(int) * (size_t) numProcs);
1921  */
1922 
1923 
1924  newCom.alltoallv->sdispls[0] = 0;
1925  newCom.alltoallv->sendcnts[0] = 0;
1926  int g,g_ant;
1927  g = g_ant = 0;
1928  for(i=0;i<nSets;i++){
1929 
1930  g = hit_lgr_elementGroup(lay,i);
1931 
1932  if(g != g_ant){
1933  newCom.alltoallv->sdispls[g] = newCom.alltoallv->sdispls[g-1] + newCom.alltoallv->sendcnts[g-1];
1934  newCom.alltoallv->sendcnts[g] = 0;
1935  g_ant = g;
1936  }
1937  newCom.alltoallv->sendcnts[g] += count[i];
1938 
1939  }
1940 
1941  while( g < numProcs ){
1942  g++;
1943  newCom.alltoallv->sendcnts[g] = 0;
1944  }
1945 
1946 
1947 
1948 
1949 
1950 
1951  int ok = MPI_Alltoall(newCom.alltoallv->sendcnts, 1, MPI_INT,newCom.alltoallv->recvcnts, 1, MPI_INT,comm);
1952  hit_mpiTestError(ok,"hit_comAllDistribute: Failed Alltoall");
1953 
1954  newCom.alltoallv->rdispls[0] = 0;
1955  for(i=1; i<numProcs; i++)
1956  newCom.alltoallv->rdispls[i] = newCom.alltoallv->rdispls[i-1] + newCom.alltoallv->recvcnts[i-1];
1957 
1958 
1959 
1960 
1961  const HitTile *tileSend = (const HitTile *)tilePSend;
1962  newCom.dataSend = tileSend->data;
1963 
1964  const HitTile *tileRecv = (const HitTile *)tilePRecv;
1965  newCom.dataRecv = tileRecv->data;
1966 
1967 
1968  /* 6. COMMUNICATION TYPE */
1969  newCom.commType = HIT_ALLTOALLV; //HIT_ALLDISTRIBUTE;
1970 
1971  /* 7. TYPE */
1972  newCom.typeSend = baseType;
1973 
1974 
1975  /* 8. RETURN */
1976  return newCom;
1977 
1978 }
1979 #endif
1980 
1981 
1982 // @note @javfres This is used in the heat benchmark to use
1983 // the same communication with different tiles. It is not pretty
1984 // and it only works for vertices, not edges.
1985 void hit_comUpdateOriginData(HitCom * com, const void * tileP){
1986  //@javfres 2015-10-05 Fix for inactive processes
1987  if(com->commType == HIT_COMTYPE_NULL) return;
1988  const HitTile *tile = (const HitTile *)tileP;
1989  com->alltoallv->sparse->originData = tile->dataVertices;
1990 }
1991 
1992 
1993 
1994 
1995 void hit_comFree(HitCom issue){
1996 
1997  switch ( issue.commType) {
1998 
1999  case HIT_SENDRECV:
2000  case HIT_REDUCE:
2001  case HIT_BROADCAST:
2002  case HIT_SENDRECV_REPLACE:
2003  hit_comFreeType(issue.typeSend);
2004  hit_comFreeType(issue.typeRecv);
2005  break;
2006 
2007  case HIT_ALLTOALL:
2008  break;
2009 
2010  case HIT_ALLTOALLV:
2011  free(issue.alltoallv->sendcnts);
2012  free(issue.alltoallv->sdispls);
2013  free(issue.alltoallv->recvcnts);
2014  free(issue.alltoallv->rdispls);
2015  free(issue.alltoallv);
2016  break;
2017  // @todo @javfres Free the lists that are now created in the communication.
2018  case HIT_SPARSE_UPDATE:
2019  free(issue.dataSend);
2020  free(issue.dataRecv);
2021  free(issue.alltoallv->sparse);
2022  free(issue.alltoallv);
2023  break;
2024  case HIT_SPARSE_SCATTER:
2025  if(hit_Rank == 0){
2026  free(issue.alltoallv->sendcnts);
2027  free(issue.alltoallv->sdispls);
2028  free(issue.dataSend);
2029  }
2030  free(issue.alltoallv->sparse);
2031  free(issue.alltoallv);
2032  break;
2033  case HIT_ALLGATHERV:
2034  free(issue.alltoallv->recvcnts);
2035  free(issue.alltoallv->rdispls);
2036  free(issue.alltoallv);
2037  break;
2038  default:
2039  break;
2040  };
2041 
2042 }
2043 
2044 
2045 
2046 
2047 
2048 
2049 
2050 /* Hit START COMMITED SEND COMMUNICATION */
2051 void hit_comStartSend(HitCom *issue) {
2052 #ifdef DEBUG
2053  printf("CTRL %d - Sends (async) from %d to %d\n", hit_Rank, issue->myself, issue->sendTo);
2054 #endif
2055 
2056  if ( issue->sendTo == MPI_PROC_NULL ) return;
2057  int ok = MPI_Isend( issue->dataSend, 1, issue->typeSend, issue->sendTo, issue->tag, issue->comm, &(issue->requestSend) );
2058  hit_mpiTestError(ok,"Failed start send");
2059 }
2060 
2061 /* Hit END COMMITED SEND COMMUNICATION */
2062 void hit_comEndSend(HitCom *issue) {
2063 #ifdef DEBUG
2064  printf("CTRL %d - End Send from %d to %d\n", hit_Rank, issue->myself, issue->sendTo);
2065 #endif
2066 
2067  if ( issue->sendTo == MPI_PROC_NULL ) return;
2068  int ok = MPI_Wait( &(issue->requestSend), MPI_STATUS_IGNORE );
2069  hit_mpiTestError(ok,"Failed end send");
2070 }
2071 
2072 /* Hit COMMITED SEND COMMUNICATION */
2073 void hit_comDoSend(HitCom *issue) {
2074 #ifdef DEBUG
2075  printf("CTRL %d - Sends from %d to %d\n", hit_Rank, issue->myself, issue->sendTo);
2076 #endif
2077  if ( issue->sendTo == MPI_PROC_NULL ) return;
2078  int ok = MPI_Send( issue->dataSend, 1, issue->typeSend, issue->sendTo, issue->tag, issue->comm );
2079  hit_mpiTestError(ok,"Failed send");
2080 }
2081 
2082 
2083 /* Hit COMMITED RECEIVE COMMUNICATION */
2084 void hit_comDoRecv(HitCom *issue) {
2085 #ifdef DEBUG
2086  printf("CTRL %d - Receives from %d to %d\n", hit_Rank, issue->recvFrom, issue->myself);
2087 #endif
2088  if ( issue->recvFrom == MPI_PROC_NULL ) return;
2089  int ok = MPI_Recv( issue->dataRecv, 1, issue->typeRecv, issue->recvFrom, issue->tag, issue->comm, MPI_STATUS_IGNORE );
2090  hit_mpiTestError(ok,"Failed receive");
2091 }
2092 
2093 /* Hit START COMMITED RECEIVE COMMUNICATION */
2094 void hit_comStartRecv(HitCom *issue) {
2095 #ifdef DEBUG
2096  printf("CTRL %d - Receives (async) from %d to %d\n", hit_Rank, issue->recvFrom, issue->myself);
2097 #endif
2098  if ( issue->recvFrom == MPI_PROC_NULL ) return;
2099  int ok = MPI_Irecv( issue->dataRecv, 1, issue->typeRecv, issue->recvFrom, issue->tag, issue->comm, &(issue->requestRecv) );
2100  hit_mpiTestError(ok,"Failed start receive");
2101 }
2102 
2103 /* Hit END COMMITED RECV COMMUNICATION */
2104 void hit_comEndRecv(HitCom *issue) {
2105 #ifdef DEBUG
2106  printf("CTRL %d - End recv from %d to %d\n", hit_Rank, issue->recvFrom, issue->myself);
2107 #endif
2108 
2109  if ( issue->sendTo == MPI_PROC_NULL ) return;
2110  int ok = MPI_Wait( &(issue->requestRecv), MPI_STATUS_IGNORE );
2111  hit_mpiTestError(ok,"Failed end recv");
2112 }
2113 
2114 
2115 /* Hit DO COMMITED SEND_RECV_REPLACE */
2117 
2118  /* 2. DO COMMUNICATION */
2119  int ok;
2120  ok = MPI_Sendrecv_replace(issue->dataSend, 1, issue->typeSend, issue->sendTo, issue->tag, issue->recvFrom ,issue->tag, issue->comm, MPI_STATUS_IGNORE);
2121  hit_mpiTestError(ok,"Failed send-receive-replace");
2122 }
2123 
2124 
2125 /* Hit DO COMMITED REDUCE */
2126 void hit_comDoReduce(HitCom *issue) {
2127 
2128 #ifdef DEBUG
2129  printf("Do reduce to %d(%d)\n",issue->sendTo,MPI_PROC_NULL);
2130 #endif
2131 
2132  int ok;
2133  // mpiReduce is mpiAllreduce if sendTo is MPI_PROC_NULL
2134  if(issue->sendTo != MPI_PROC_NULL){
2135  ok = MPI_Reduce( issue->dataSend, issue->dataRecv, 1, issue->typeSend, issue->operation, issue->sendTo, issue->comm);
2136  } else {
2137  ok = MPI_Allreduce(issue->dataSend, issue->dataRecv, 1, issue->typeSend, issue->operation, issue->comm);
2138  }
2139  hit_mpiTestError(ok,"Failed reduction");
2140 }
2141 
2142 /* Hit DO COMMITED BROADCAST */
2144 
2145 #ifdef DEBUG
2146  printf("CTRL %d - Broadcast from %d\n", hit_Rank, issue->sendTo);
2147 #endif
2148 
2149  int ok;
2150  ok = MPI_Bcast( issue->dataSend, 1, issue->typeSend, issue->sendTo, issue->comm);
2151  hit_mpiTestError(ok,"Failed broadcast");
2152 }
2153 
2154 /* Hit DO COMMITED ALLTOALL */
2156 
2157  int ok;
2158  ok = MPI_Alltoall(issue->dataSend, issue->count, issue->typeSend, issue->dataRecv, issue->count, issue->typeSend, issue->comm);
2159 
2160  hit_mpiTestError(ok,"Failed all-to-all");
2161 }
2162 
2163 
2164 /* Hit DO COMMITED Alltoallv */
2166 
2167  int ok;
2168  ok = MPI_Alltoallv(issue->dataSend, issue->alltoallv->sendcnts, issue->alltoallv->sdispls,issue->typeSend, issue->dataRecv, issue->alltoallv->recvcnts,issue->alltoallv->rdispls, issue->typeSend, issue->comm);
2169 
2170  hit_mpiTestError(ok,"Failed all-to-all v");
2171 }
2172 
2173 /* Hit DO COMMITED SparseUpdate for CSR shapes */
2175 
2176 
2177  int i;
2178  size_t baseExtent = issue->alltoallv->sparse->baseExtent;
2179 
2180 
2181  /* 1. Copy the local values to the send vector */
2182  for(i=0;i<issue->alltoallv->sparse->nsend;i++){
2183  int vertex = issue->alltoallv->sparse->send[i];
2184  int lvertex = hit_cShapeVertexToLocal(issue->alltoallv->sparse->shape,vertex);
2185 
2186  memcpy( (char*)issue->dataSend + (size_t) i * baseExtent,
2187  (char*)issue->alltoallv->sparse->originData + (size_t) lvertex * baseExtent ,
2188  (size_t)baseExtent);
2189 
2190  //debugall("Sending %d->: %d, local %d, value %f\n",issue->myself,vertex,lvertex, *((double*)(issue->dataSend + i * baseExtent)));
2191  }
2192 
2193 
2194  /* 2. All to all communication. */
2195  int ok = MPI_Alltoallv(issue->dataSend, issue->alltoallv->sendcnts,
2196  issue->alltoallv->sdispls,issue->typeSend,
2197  issue->dataRecv, issue->alltoallv->recvcnts,issue->alltoallv->rdispls,
2198  issue->typeSend, issue->comm);
2199 
2200  hit_mpiTestError(ok,"Failed sparse update");
2201 
2202  /* 3. Extract the updated values from the recv vector. */
2203  for(i=0;i<issue->alltoallv->sparse->nrecv;i++){
2204  int vertex = issue->alltoallv->sparse->recv[i];
2205 
2206  int lvertex = hit_cShapeVertexToLocal(issue->alltoallv->sparse->shape,vertex);
2207 
2208  memcpy( (char*)issue->alltoallv->sparse->originData + (size_t) lvertex * baseExtent,
2209  (char*)issue->dataRecv + (size_t) i * baseExtent,
2210  baseExtent);
2211 
2212  //debugall("Recving %d<-: %d, local %d, value %d\n",issue->myself,vertex,lvertex, *((int*)(issue->dataRecv + i * baseExtent)));
2213  }
2214 }
2215 
2216 /* Hit DO COMMITED SparseUpdate for Bitmap shapes */
2218 
2219  int i;
2220  size_t baseExtent = issue->alltoallv->sparse->baseExtent;
2221 
2222 
2223  /* 1. Copy the local values to the send vector */
2224  for(i=0;i<issue->alltoallv->sparse->nsend;i++){
2225  int vertex = issue->alltoallv->sparse->send[i];
2226  int lvertex = hit_bShapeVertexToLocal(issue->alltoallv->sparse->shape,vertex);
2227 
2228  memcpy( (char*)issue->dataSend + (size_t) i * baseExtent,
2229  (char*)issue->alltoallv->sparse->originData + (size_t) lvertex * baseExtent ,
2230  (size_t)baseExtent);
2231 
2232  //printf("Sending %d->: %d, local %d, value %d\n",issue->myself,vertex,lvertex, *((int*)(issue->dataSend + i * baseExtent)));
2233  }
2234 
2235 
2236  /* 2. All to all communication. */
2237  int ok = MPI_Alltoallv(issue->dataSend, issue->alltoallv->sendcnts,
2238  issue->alltoallv->sdispls,issue->typeSend,
2239  issue->dataRecv, issue->alltoallv->recvcnts,issue->alltoallv->rdispls,
2240  issue->typeSend, issue->comm);
2241 
2242  hit_mpiTestError(ok,"Failed sparse update");
2243 
2244  /* 3. Extract the updated values from the recv vector. */
2245  for(i=0;i<issue->alltoallv->sparse->nrecv;i++){
2246  int vertex = issue->alltoallv->sparse->recv[i];
2247 
2248  int lvertex = hit_bShapeVertexToLocal(issue->alltoallv->sparse->shape,vertex);
2249 
2250  memcpy( (char*)issue->alltoallv->sparse->originData + (size_t) lvertex * baseExtent,
2251  (char*)issue->dataRecv + (size_t) i * baseExtent,
2252  baseExtent);
2253 
2254  //debugall("Recving %d<-: %d, local %d, value %d\n",issue->myself,vertex,lvertex, *((int*)(issue->dataRecv + i * baseExtent)));
2255  }
2256 }
2257 
2258 /* Hit DO COMMITED SparseUpdate */
2260 
2261  int type = hit_shapeType(issue->alltoallv->sparse->shape);
2262 
2263  switch(type){
2264  case HIT_CSR_SHAPE:
2265  hit_comDoSparseUpdateCSR(issue);
2266  break;
2267  case HIT_BITMAP_SHAPE:
2269  break;
2270  default:
2271  hit_errInternal(__func__, "Unsupported shape type", "", __FILE__, __LINE__);
2272  break;
2273  }
2274 }
2275 
2277 
2278  int mpi_rank;
2279  MPI_Comm_rank(issue->comm,&mpi_rank);
2280 
2281  if(mpi_rank == 0){
2282 
2283  size_t baseExtent = issue->alltoallv->sparse->baseExtent;
2284 
2285  int nProcs;
2286  MPI_Comm_size(issue->comm,&nProcs);
2287  int nElems = issue->alltoallv->sparse->nsend;
2288  int * assignedGroups = issue->alltoallv->sparse->send;
2289 
2290  int * count_aux;
2291  // @arturo Ago 2015: New allocP interface
2292  // hit_calloc(count_aux,(size_t) nProcs, sizeof(int),int*);
2293  hit_calloc(count_aux, int, nProcs);
2294 
2295  int i;
2296  for(i=0;i<nElems;i++){
2297 
2298  int dst = assignedGroups[i];
2299 
2300  memcpy( (char*)issue->dataSend + (size_t) ( issue->alltoallv->sdispls[dst] + count_aux[dst]) * baseExtent,
2301  (char*)issue->alltoallv->sparse->originData + (size_t) i * baseExtent ,
2302  (size_t)baseExtent);
2303 
2304 
2305  count_aux[ dst ]++;
2306  }
2307 
2308 
2309  free(count_aux);
2310 
2311  } else {
2312 
2313  issue->dataSend = issue->alltoallv->sendcnts = issue->alltoallv->sdispls = NULL;
2314 
2315  }
2316 
2317  int ok = MPI_Scatterv( issue->dataSend, issue->alltoallv->sendcnts, issue->alltoallv->sdispls,
2318  issue->typeSend , issue->dataRecv , issue->count,
2319  issue->typeRecv,
2320  0, issue->comm);
2321  hit_mpiTestError(ok,"Failed sparse scatter");
2322 
2323 }
2324 
2325 
2326 
2328 
2329  int mpi_rank;
2330  MPI_Comm_rank(issue->comm,&mpi_rank);
2331 
2332  if(mpi_rank == 0){
2333 
2334  size_t baseExtent = issue->alltoallv->sparse->baseExtent;
2335 
2336  int nProcs;
2337  MPI_Comm_size(issue->comm,&nProcs);
2338  int nRows = issue->alltoallv->sparse->nsend;
2339  int * assignedGroups = issue->alltoallv->sparse->send;
2340 
2341  int * count_aux;
2342  // @arturo Ago 2015: New allocP interface
2343  // hit_calloc(count_aux,(size_t) nProcs, sizeof(int),int*);
2344  hit_calloc(count_aux, int, nProcs);
2345 
2346  int i;
2347  for(i=0;i<nRows;i++){
2348 
2349  int dst = assignedGroups[i];
2350  int begin = issue->alltoallv->sparse->rows[i];
2351  int end = issue->alltoallv->sparse->rows[i+1];
2352 
2353  int j;
2354  for(j=begin; j<end; j++){
2355 
2356  memcpy( (char*)issue->dataSend + (size_t) ( issue->alltoallv->sdispls[dst] + count_aux[dst]) * baseExtent,
2357  (char*)issue->alltoallv->sparse->originData + (size_t) j * baseExtent ,
2358  (size_t)baseExtent);
2359  count_aux[ dst ]++;
2360 
2361  }
2362 
2363 
2364  //memcpy( (char*)issue->dataSend + (size_t) ( issue->alltoallv->sdispls[dst] + count_aux[dst]) * baseExtent,
2365  // (char*)issue->alltoallv->sparse->originData + (size_t) i * baseExtent ,
2366  // (size_t)baseExtent);
2367 
2368 
2369  //count_aux[ dst ]++;
2370  }
2371 
2372 
2373  free(count_aux);
2374 
2375  } else {
2376 
2377  issue->dataSend = issue->alltoallv->sendcnts = issue->alltoallv->sdispls = NULL;
2378 
2379  }
2380 
2381  int ok = MPI_Scatterv( issue->dataSend, issue->alltoallv->sendcnts, issue->alltoallv->sdispls,
2382  issue->typeSend , issue->dataRecv , issue->count,
2383  issue->typeRecv,
2384  0, issue->comm);
2385  hit_mpiTestError(ok,"Failed sparse scatter");
2386 
2387 }
2388 
2389 
2390 
2392 
2393  int ok;
2394  ok = MPI_Allgatherv(issue->dataSend,issue->count,issue->typeSend,
2395  issue->dataRecv,
2396  issue->alltoallv->recvcnts,
2397  issue->alltoallv->rdispls,
2398  issue->typeRecv,issue->comm);
2399  hit_mpiTestError(ok,"Failed all_gatherv");
2400 }
2401 
2402 
2403 
2404 
2405 
2406 
2407 /* hit_comDo: executes a Hit Communication issue */
2408 void hit_comDo(HitCom *issue) {
2409  /* 1. SKIP NULL COMM */
2410  if ( issue->commType == HIT_COMTYPE_NULL ) return;
2411 
2412  /* 2. CALL THE PROPER FUNCTION TO DO THE PATTERN */
2413  switch ( issue->commType) {
2414  case HIT_SENDRECV:
2415  /* a. START COUPLED ISEND */
2416  hit_comStartSend(issue);
2417  /* b. ACCEPT COUPLED RECEIVE */
2418  hit_comDoRecv(issue);
2419  /* c. WAIT FOR COUPLED ISEND TO COMPLETE */
2420  hit_comEndSend(issue);
2421  break;
2422  case HIT_REDUCE:
2423  hit_comDoReduce(issue);
2424  break;
2425  case HIT_ALLTOALL:
2426  hit_comDoAlltoall(issue);
2427  break;
2428  case HIT_ALLTOALLV:
2429  hit_comDoAlltoallv(issue);
2430  break;
2431  case HIT_BROADCAST:
2432  hit_comDoBroadcast(issue);
2433  break;
2434  case HIT_SENDRECV_REPLACE:
2435  hit_comDoSendRecvReplace(issue);
2436  break;
2437  case HIT_SPARSE_UPDATE:
2438  hit_comDoSparseUpdate(issue);
2439  break;
2440  case HIT_SPARSE_SCATTER:
2441  hit_comDoSparseScatter(issue);
2442  break;
2445  break;
2446  case HIT_ALLGATHERV:
2447  hit_comDoAllGatherv(issue);
2448  break;
2449  default:
2450  hit_errInternal("comDo", "Unknown type of communication object", "", __FILE__, __LINE__);
2451  break;
2452  };
2453 }
2454 
2460 void hit_comOpSumDoubleBasic(void * in, void * inout) {
2461 
2462  double * din = (double *)in;
2463  double * dinout = (double *)inout;
2464 
2465 
2466 #ifdef DEBUG
2467  printf("CTRL OpSumDouble: %lf+=%lf=%lf \n",*dinout,*din,*dinout+*din);
2468 #endif
2469 
2470  *dinout+=*din;
2471 }
2472 
2478 void hit_comOpMaxDoubleBasic(void * in, void * inout) {
2479 
2480  double * din = (double *)in;
2481  double * dinout = (double *)inout;
2482 
2483 #ifdef DEBUG
2484  printf("CTRL OpMaxDouble: hit_max(%lf,%lf)=%lf \n",*dinout,*din,hit_max(*dinout,*din));
2485 #endif
2486 
2487  *dinout=hit_max(*dinout,*din);
2488 }
2489 
2495 void hit_comOpMinDoubleBasic(void * in, void * inout) {
2496 
2497  double * din = (double *)in;
2498  double * dinout = (double *)inout;
2499 
2500 #ifdef DEBUG
2501  printf("CTRL OpMinDouble: hit_min(%lf,%lf)=%lf \n",*dinout,*din,hit_min(*dinout,*din));
2502 #endif
2503 
2504  *dinout=hit_min(*dinout,*din);
2505 }
2506 
2512 void hit_comOpSumIntBasic(void * in, void * inout) {
2513 
2514  int * din = (int *)in;
2515  int * dinout = (int *)inout;
2516 
2517  *dinout+=*din;
2518 }
2519 
2520 
2526 void hit_comOpMaxIntBasic(void * in, void * inout) {
2527 
2528  int * din = (int *)in;
2529  int * dinout = (int *)inout;
2530 
2531  *dinout=hit_max(*dinout,*din);
2532 }
2533 
2534 
2535 
2541 void hit_comOpMinIntBasic(void * in, void * inout) {
2542 
2543  int * din = (int *)in;
2544  int * dinout = (int *)inout;
2545 
2546  *dinout=hit_min(*dinout,*din);
2547 }
2548 
2549 
2550 
2551 /* generic function to make operations between MPI user datatypes */
2552 void hit_comOpGenericAnyType( void * in, void * inout, MPI_Datatype datatype, int offset, size_t tam, HitComOpFunction f )
2553 {
2554  int *array_of_ints;
2555  MPI_Aint *array_of_adds;
2556  MPI_Datatype *array_of_dtypes;
2557  int num_ints, num_adds, num_dtypes, combiner;
2558  int i,j,stride,jump,ac,ac2;
2559  MPI_Aint ext;
2560 
2561  MPI_Type_get_envelope( datatype,
2562  &num_ints, &num_adds, &num_dtypes, &combiner );
2563 
2564  switch (combiner) {
2565  case MPI_COMBINER_NAMED: // Native MPI type
2566  stride=offset/(int)tam;
2567  f( (char *)in+offset, (char *)inout+offset);
2568  break;
2569  case MPI_COMBINER_STRUCT: //N+1-N-N N is number of datatypes, stored in ints[0]
2570 
2571  hit_malloc(array_of_ints, int, num_ints);
2572  hit_malloc(array_of_adds, MPI_Aint, num_adds);
2573  hit_malloc(array_of_dtypes, MPI_Datatype, num_dtypes);
2574  // @arturo Ago 2015: New allocP interface
2575  /*
2576  hit_malloc(array_of_ints,(size_t) num_ints * sizeof(int),int*);
2577  hit_malloc(array_of_adds,(size_t) num_adds * sizeof(MPI_Aint),MPI_Aint*);
2578  hit_malloc(array_of_dtypes,(size_t) num_dtypes * sizeof(MPI_Datatype),MPI_Datatype*);
2579  */
2580 
2581 
2582  MPI_Type_get_contents( datatype, num_ints, num_adds, num_dtypes,
2583  array_of_ints, array_of_adds, array_of_dtypes );
2584  //our blocklengths are always one: internal loop not tested
2585  for (i=0; i<array_of_ints[0]; i++) {
2586  ac= (int) array_of_adds[i] + offset;
2587  MPI_Type_extent(array_of_dtypes[i],&ext);
2588  jump = (int) ext;
2589  for(j=0;j<array_of_ints[i+1];j++) {
2590  hit_comOpGenericAnyType( in, inout, array_of_dtypes[i],ac,tam,f);
2591  ac+=jump;
2592  }
2593  }
2594  free( array_of_ints );
2595  free( array_of_adds );
2596  free( array_of_dtypes );
2597  break;
2598  case MPI_COMBINER_HVECTOR: //2-1-1 count in ints[0] & blocklength in ints[1] & stride in adds[0] & oldtype in dtypes[0]
2599 
2600 
2601  hit_malloc(array_of_ints, int, num_ints);
2602  hit_malloc(array_of_adds, MPI_Aint, num_adds);
2603  hit_malloc(array_of_dtypes, MPI_Datatype, num_dtypes);
2604  // @arturo Ago 2015: New allocP interface
2605  /*
2606  hit_malloc(array_of_ints,(size_t) num_ints * sizeof(int), int* );
2607  hit_malloc(array_of_adds,(size_t) num_adds * sizeof(MPI_Aint), MPI_Aint*);
2608  hit_malloc(array_of_dtypes,(size_t) num_dtypes * sizeof(MPI_Datatype),MPI_Datatype*);
2609  */
2610 
2611 
2612  MPI_Type_get_contents( datatype, num_ints, num_adds, num_dtypes,
2613  array_of_ints, array_of_adds, array_of_dtypes );
2614  //our blocklengths are always one: internal loop not tested
2615  MPI_Type_extent(array_of_dtypes[0],&ext);
2616  jump = (int) ext;
2617  // @javfres ac2 does nothing
2618  //ac=offset, ac2;
2619  ac=offset;
2620  stride = (int) array_of_adds[0];
2621  for(i=0;i<array_of_ints[0];i++) {
2622  ac2=ac;
2623  for(j=0;j<array_of_ints[1];j++) {
2624  hit_comOpGenericAnyType( in, inout, array_of_dtypes[0],ac2, tam,f );
2625  ac2+=jump;
2626  }
2627  ac+=stride;
2628  }
2629  free( array_of_ints );
2630  free( array_of_adds );
2631  free( array_of_dtypes );
2632  break;
2633  case MPI_COMBINER_CONTIGUOUS: //1-0-1 count in ints[0] & oldtype in dtypes[0]
2634 
2635  hit_malloc(array_of_ints, int, num_ints);
2636  hit_malloc(array_of_adds, MPI_Aint, num_adds);
2637  hit_malloc(array_of_dtypes, MPI_Datatype, num_dtypes);
2638  // @arturo Ago 2015: New allocP interface
2639  /*
2640  hit_malloc(array_of_ints,(size_t) num_ints * sizeof(int) ,int* );
2641  hit_malloc(array_of_adds,(size_t) num_adds * sizeof(MPI_Aint),MPI_Aint*);
2642  hit_malloc(array_of_dtypes,(size_t) num_dtypes * sizeof(MPI_Datatype),MPI_Datatype*);
2643  */
2644 
2645 
2646  MPI_Type_get_contents( datatype, num_ints, num_adds, num_dtypes,
2647  array_of_ints, array_of_adds, array_of_dtypes );
2648  MPI_Type_extent(array_of_dtypes[0],&ext);
2649  jump = (int) ext;
2650  ac=offset;
2651  for(i=0;i<array_of_ints[0];i++) {
2652  hit_comOpGenericAnyType( in, inout, array_of_dtypes[0],ac,tam,f );
2653  ac+=jump;
2654  }
2655  free( array_of_ints );
2656  free( array_of_adds );
2657  free( array_of_dtypes );
2658  break;
2659  }
2660 }
2661 
2662 
2663 /* Hit OPERATION to add doubles in reductions */
2664 //@see decoding a datatype in mpi_forum.org
2665 void hit_comOpSumDouble (void * in, void * inout, int *len, HitType * type) {
2666  HIT_NOT_USED(len);
2667  hit_comOpGenericAnyType(in,inout,*type,0,sizeof(double),&hit_comOpSumDoubleBasic);
2668 }
2669 
2670 /* Hit OPERATION to calculate the maximum in reductions */
2671 void hit_comOpMaxDouble (void * in, void * inout, int *len, HitType * type) {
2672  HIT_NOT_USED(len);
2673  hit_comOpGenericAnyType(in,inout,*type,0,sizeof(double),&hit_comOpMaxDoubleBasic);
2674 }
2675 
2676 /* Hit OPERATION to calculate the minimum in reductions */
2677 void hit_comOpMinDouble (void * in, void * inout, int *len, HitType * type) {
2678  HIT_NOT_USED(len);
2679  hit_comOpGenericAnyType(in,inout,*type,0,sizeof(double),&hit_comOpMinDoubleBasic);
2680 }
2681 
2682 
2683 void hit_comOpSumInt (void * in, void * inout, int *len, HitType * type) {
2684  HIT_NOT_USED(len);
2685  hit_comOpGenericAnyType(in,inout,*type,0,sizeof(int),&hit_comOpSumIntBasic);
2686 }
2687 void hit_comOpMaxInt (void * in, void * inout, int *len, HitType * type) {
2688  HIT_NOT_USED(len);
2689  hit_comOpGenericAnyType(in,inout,*type,0,sizeof(int),&hit_comOpMaxIntBasic);
2690 }
2691 void hit_comOpMinInt (void * in, void * inout, int *len, HitType * type) {
2692  HIT_NOT_USED(len);
2693  hit_comOpGenericAnyType(in,inout,*type,0,sizeof(int),&hit_comOpMinIntBasic);
2694 }
2695 
2696 
2697 
void hit_comDoSend(HitCom *issue)
Definition: hit_com.c:2073
HitCom hit_comSparseUpdateBitmap(HitLayout lay, const void *tileP, HitType baseType)
Definition: hit_com.c:1324
HitOp HIT_OP_MIN_DOUBLE
Definition: hit_com.c:67
void hit_comDoReduce(HitCom *issue)
Definition: hit_com.c:2126
MPI_Comm comm
Definition: hit_com.h:235
Vector g(int i, int j)
MPI_Op HitOp
Definition: hit_com.h:118
int * sizes
Definition: refMPIluBack.c:104
HitType hit_comTypeRec(const void *varP, HitType baseType)
Definition: hit_com.c:290
struct HitComSparse * sparse
Definition: hit_com.h:259
#define hit_free(ptr)
Definition: hit_allocP.h:152
#define hit_layShape(lay)
Definition: hit_layout.h:650
int * send
Definition: hit_com.h:273
#define hit_tileDimStride(var, dim)
Definition: hit_tile.h:806
HitCom hit_comSparseUpdate(HitLayout lay, const void *tileP, HitType baseType)
Definition: hit_com.c:1084
HitComAlltoallv HIT_COM_ALLTOALLV_NULL
Definition: hit_com.c:59
void hit_comOpSumDouble(void *, void *, int *, HitType *)
Definition: hit_com.c:2665
void hit_comOpMinIntBasic(void *in, void *inout)
Definition: hit_com.c:2541
MPI_Request requestSend
Definition: hit_com.h:236
#define hit_comOpFree(operation)
Definition: hit_com.h:1048
#define HIT_ALLGATHERV
Definition: hit_com.h:113
HitOp HIT_OP_MAX_DOUBLE
Definition: hit_com.c:68
HitComSparse HIT_COM_SPARSE_NULL
Definition: hit_com.c:60
HitOp HIT_OP_SUM_INT
Definition: hit_com.c:63
HitPTopology * hit_ptopSplit(HitPTopology *in, int group)
Definition: hit_topology.c:80
HitCom hit_comBroadcastSelect(HitLayout lay, HitRanks root, const void *tile, HitShape selection, int modeSelect, HitType baseType)
Definition: hit_com.c:731
#define hit_Rank
Definition: hit_com.h:140
HitCom hit_comAllGathervInternal(HitLayout lay, const void *tilePSend, const void *tilePRecv, HitType baseType, const char *file, int line)
Definition: hit_com.c:996
int hit_topRankInternal(HitTopology topo, HitRanks ranks)
Definition: hit_topology.c:415
void hit_comOpMaxDoubleBasic(void *in, void *inout)
Definition: hit_com.c:2478
void hit_comOpMaxInt(void *, void *, int *, HitType *)
Definition: hit_com.c:2687
void hit_comDoSendRecvReplace(HitCom *issue)
Definition: hit_com.c:2116
int rank[HIT_MAXDIMS]
Definition: hit_topology.h:135
HitTopology topo
Definition: hit_layout.h:255
void hit_comEndRecv(HitCom *issue)
Definition: hit_com.c:2104
#define hit_laySelfRanks(lay)
Definition: hit_layout.h:848
HitPTopology * pTopology
Definition: hit_topology.h:253
HitOp HIT_OP_MIN_INT
Definition: hit_com.c:64
HitPTopology * HIT_TOPOLOGY_INFO
Definition: hit_topology.c:60
#define hit_tileDims(var)
Definition: hit_tile.h:713
void hit_comFree(HitCom issue)
Definition: hit_com.c:1995
HitCom hit_comSendRecvSelectTag(HitLayout lay, HitRanks sendTo, const void *tilePSend, HitShape selectionSend, int modeSelectSend, HitRanks receiveFrom, const void *tilePRecv, HitShape selectionRecv, int modeSelectRecv, HitType baseType, int tag)
Definition: hit_com.c:450
int numElementsTotal
Definition: hit_layout.h:226
void hit_comDo(HitCom *issue)
Definition: hit_com.c:2408
int count
Definition: hit_com.h:242
HitType HIT_SHAPE_SIG
Definition: hit_com.c:71
#define hit_bShapeNameList(shape, dim)
Definition: hit_bshape.h:171
void hit_comStartRecv(HitCom *issue)
Definition: hit_com.c:2094
#define hit_tileDimSig(var, dim)
Definition: hit_tile.h:776
#define hit_cShapeNColsRow(s, row)
Definition: hit_cshape.h:210
#define HIT_COM_SPARSE_NULL_STATIC
Definition: hit_com.h:335
void hit_comDoSparseScatter(HitCom *issue)
Definition: hit_com.c:2276
HitCom hit_comSparseScatterInternal(HitLayout lay, const void *tilePSend, const void *tilePRecv, HitType baseType, const char *file, int line)
Definition: hit_com.c:1543
void hit_comDoSparseScatterRows(HitCom *issue)
Definition: hit_com.c:2327
#define hit_calloc(ptr, type, nmemb)
Definition: hit_allocP.h:114
void hit_comOpMaxIntBasic(void *in, void *inout)
Definition: hit_com.c:2526
int hit_tileCheckBoundaryArrayCoords(const void *tileP, HitShape sh)
Definition: hit_tile.c:1073
#define HIT_SPARSE_UPDATE
Definition: hit_com.h:105
HitCom hit_comReduceSelect(HitLayout lay, HitRanks root, const void *tilePSend, HitShape selectionSend, int modeSelectSend, const void *tilePRecv, HitShape selectionRecv, int modeSelectRecv, HitType baseType, HitOp operation)
Definition: hit_com.c:585
#define HIT_ALLTOALL
Definition: hit_com.h:97
HitType typeSend
Definition: hit_com.h:231
void hit_comDoSparseUpdateBitmap(HitCom *issue)
Definition: hit_com.c:2217
HitCom hit_comAllDistribute(HitLayout lay, const void *tilePSend, int *count, const void *tilePRecv, HitType baseType)
Definition: hit_com.c:1651
HitCom hit_comSparseUpdateCSR(HitLayout lay, const void *tileP, HitType baseType)
Definition: hit_com.c:1109
void hit_comOpSumIntBasic(void *in, void *inout)
Definition: hit_com.c:2512
#define hit_layShapeOther(lay, ranks)
Definition: hit_layout.h:732
MPI_Request requestRecv
Definition: hit_com.h:237
HitRanks HIT_RANKS_NULL
Definition: hit_topology.c:68
#define hit_sigCard(sig)
Definition: hit_sig.h:162
int * rows
Definition: hit_com.h:279
#define hit_ranksCmp(a, b)
Definition: hit_topology.h:220
#define HIT_MAXDIMS
Definition: hit_shape.h:72
void * dataSend
Definition: hit_com.h:233
void hit_comEndSend(HitCom *issue)
Definition: hit_com.c:2062
int HIT_RANK_NULL
Definition: hit_topology.c:67
HitPTopology * pTopology[HIT_MAXDIMS+1]
Definition: hit_layout.h:278
int hit_layActiveRanksId(HitLayout lay, HitRanks ranks)
Definition: hit_layout.c:1916
HitCom hit_comAlltoallSelectv(HitLayout lay, const void *tilePSend, HitShape *selectionSend, int modeSelectSend, const void *tilePRecv, HitShape *selectionRecv, int modeSelectRecv, HitType baseType)
Definition: hit_com.c:896
void * originData
Definition: hit_com.h:270
#define hit_tileCard(var)
Definition: hit_tile.h:763
int * sendcnts
Definition: hit_com.h:254
int nProcs
Definition: SWpar_ref.c:183
#define HIT_COM_ARRAYCOORDS
Definition: hit_com.h:345
#define HIT_SPARSE_SCATTER
Definition: hit_com.h:107
#define HIT_BITMAP_SHAPE
Definition: hit_shape.h:183
MPI_Datatype * types
Definition: refMPIluBack.c:109
void hit_comOpSumInt(void *, void *, int *, HitType *)
Definition: hit_com.c:2683
#define hit_tileSelect(newVar, oldVar, shape)
Definition: hit_tile.h:453
int * recvcnts
Definition: hit_com.h:256
int sendTo
Definition: hit_com.h:228
char * get_gdb_trace()
Definition: hit_error.c:45
#define hit_layNumDims(lay)
Definition: hit_layout.h:655
MPI_Comm comm
Definition: SWpar_ref.c:193
void(* HitComOpFunction)(void *, void *)
Definition: hit_com.h:1092
void hit_comDoAlltoallv(HitCom *issue)
Definition: hit_com.c:2165
void hit_comDoAlltoall(HitCom *issue)
Definition: hit_com.c:2155
#define hit_tileDimCard(var, dim)
Definition: hit_tile.h:750
int * recv
Definition: hit_com.h:272
HitCom hit_comBroadcastDimSelect(HitLayout lay, int dim, int root, const void *tile, HitShape selection, int modeSelect, HitType baseType)
Definition: hit_com.c:779
void hit_comUpdateOriginData(HitCom *com, const void *tileP)
Definition: hit_com.c:1985
#define hit_shapeDims(shape)
Definition: hit_sshape.h:364
int * sdispls
Definition: hit_com.h:255
#define hit_cShapeXadj(shape)
HitType hit_comType(const void *varP, HitType baseType)
Definition: hit_com.c:189
#define hit_mpiTestError(ok, cad)
Definition: hit_com.h:62
#define HIT_ALLTOALLV
Definition: hit_com.h:99
int commType
Definition: hit_com.h:226
#define HIT_COMTYPE_NULL
Definition: hit_com.h:89
#define hit_cShapeNvertices(shape)
int * assignedGroups
Definition: hit_layout.h:227
int recvFrom
Definition: hit_com.h:229
#define hit_bShapeVertexToLocal(s, vertex)
Definition: hit_bshape.h:355
int * rdispls
Definition: hit_com.h:257
#define hit_cShapeVertexToLocal(s, vertex)
Definition: hit_cshape.h:243
#define hit_error(name, file, numLine)
Definition: hit_com.h:1012
#define hit_bShapeVertexToGlobal(s, vertex)
int tag
Definition: hit_com.h:230
#define hit_comOp(function, operation)
Definition: hit_com.h:1036
#define hit_cShapeNameList(shape, dim)
Definition: hit_cshape.h:162
MPI_Datatype HitType
Definition: hit_com.h:70
#define HIT_TYPE_NULL
Definition: hit_com.h:181
#define hit_malloc(ptr, type, nmemb)
Definition: hit_allocP.h:93
#define HIT_NOT_USED(x)
Definition: hit_error.h:50
HitShape ext_shape
#define HIT_SPARSE_SCATTER_ROWS
Definition: hit_com.h:111
#define hit_cShapeNZElems(shape)
Definition: hit_cshape.h:126
void hit_comOpMinInt(void *, void *, int *, HitType *)
Definition: hit_com.c:2691
struct HitComAlltoallv * alltoallv
Definition: hit_com.h:243
void hit_ptopFree(HitPTopology **in)
Definition: hit_topology.c:109
#define HIT_CSR_SHAPE
Definition: hit_shape.h:177
#define HIT_SENDRECV
Definition: hit_com.h:91
HitCom hit_comAlltoallSelect(HitLayout lay, const void *tilePSend, HitShape selectionSend, int modeSelectSend, const void *tilePRecv, HitShape selectionRecv, int modeSelectRecv, HitType baseType, int count)
Definition: hit_com.c:834
void hit_comDoSparseUpdateCSR(HitCom *issue)
Definition: hit_com.c:2174
int hit_tileCheckBoundaryTileCoords(const void *tileP, HitShape sh)
Definition: hit_tile.c:1095
void hit_comDoAllGatherv(HitCom *issue)
Definition: hit_com.c:2391
HitShape shape
Definition: hit_com.h:269
HitOp HIT_OP_SUM_DOUBLE
Definition: hit_com.c:66
void hit_comInit(int *pargc, char **pargv[])
Definition: hit_com.c:111
HitShape shape
void * hit_comSearchData(const void *varP)
Definition: hit_com.c:278
HitLayoutList layoutList
Definition: hit_layout.h:285
#define hit_bShapeNvertices(shape)
void hit_comDoBroadcast(HitCom *issue)
Definition: hit_com.c:2143
void hit_comOpMaxDouble(void *, void *, int *, HitType *)
Definition: hit_com.c:2671
#define HIT_LAYOUT_CONTIGUOUS
Definition: hit_layout.h:381
size_t baseExtent
Definition: hit_com.h:277
#define hit_vmalloc(ptr, size)
Definition: hit_allocP.h:72
void * dataRecv
Definition: hit_com.h:234
#define HIT_COM_MYSELF
Definition: hit_com.h:135
HitCom HIT_COM_NULL
Definition: hit_com.c:58
HitShape shape
Definition: hit_layout.h:271
MPI_Aint * addresses
Definition: refMPIluBack.c:105
int myself
Definition: hit_com.h:227
HitCom hit_comSendRecvReplaceSelectTag(HitLayout lay, HitRanks sendTo, const void *tileP, HitShape selection, int modeSelect, HitRanks receiveFrom, HitType baseType, int tag)
Definition: hit_com.c:393
#define hit_shapeType(s)
Definition: hit_shape.h:266
union HitLayout::@0 info
void hit_comOpSumDoubleBasic(void *in, void *inout)
Definition: hit_com.c:2460
#define hit_tileSelectArrayCoords(newVar, oldVar, shape)
Definition: hit_tile.h:494
#define hit_lgr_elementGroup(lay, element)
Definition: hit_layout.h:938
void hit_comOpGenericAnyType(void *in, void *inout, HitType datatype, int offset, size_t tam, HitComOpFunction f)
Definition: hit_com.c:2552
HitCom hit_comSparseScatterRowsInternal(HitLayout lay, const void *tilePSend, const void *tilePRecv, HitType baseType, const char *file, int line)
Definition: hit_com.c:1764
HitShape origShape
Definition: hit_layout.h:274
HitType typeRecv
Definition: hit_com.h:232
#define HIT_SENDRECV_REPLACE
Definition: hit_com.h:103
#define HIT_REDUCE
Definition: hit_com.h:93
#define HIT_COM_NULL_STATIC
Definition: hit_com.h:308
#define hit_cShapeVertexToGlobal(s, vertex)
void hit_comDoSparseUpdate(HitCom *issue)
Definition: hit_com.c:2259
HitCom hit_comReduceDimSelect(HitLayout lay, int dim, HitRanks root, const void *tilePSend, HitShape selectionSend, int modeSelectSend, const void *tilePRecv, HitShape selectionRecv, int modeSelectRecv, HitType baseType, HitOp operation)
Definition: hit_com.c:660
#define hit_max(a, b)
Definition: hit_funcop.h:59
int active
Definition: hit_layout.h:263
#define hit_topSelfRankInternal(topo)
Definition: hit_topology.h:439
#define hit_warnInternal(routine, text, extraParam, file, numLine)
Definition: hit_error.h:69
void hit_comAllowDim(HitLayout *lay, int dim)
Definition: hit_com.c:367
int numActives[HIT_MAXDIMS]
Definition: hit_layout.h:257
HitLayout lay
Definition: heat.c:107
#define hit_shapeSig(shape, dim)
Definition: hit_sshape.h:400
HitOp operation
Definition: hit_com.h:240
#define HIT_COM_TILECOORDS
Definition: hit_com.h:343
MPI_Errhandler mpi_errhandler
Definition: hit_com.c:74
#define hit_min(a, b)
Definition: hit_funcop.h:65
#define tag(dir)
Definition: mg.c:191
#define hit_errInternal(routine, text, extraParam, file, numLine)
Definition: hit_error.h:63
#define max(a, b)
Definition: cannonAsync.c:47
#define HIT_COM_ALLTOALLV_NULL_STATIC
Definition: hit_com.h:330
HitOp HIT_OP_MAX_INT
Definition: hit_com.c:65
#define hit_sigTileToArray(sig, ind)
Definition: hit_sig.h:198
void hit_comOpMinDoubleBasic(void *in, void *inout)
Definition: hit_com.c:2495
#define hit_comFreeType(type)
Definition: hit_com.h:187
void hit_comFinalize()
Definition: hit_com.c:159
void hit_comOpMinDouble(void *, void *, int *, HitType *)
Definition: hit_com.c:2677
#define HIT_BROADCAST
Definition: hit_com.h:101
void hit_comDoRecv(HitCom *issue)
Definition: hit_com.c:2084
#define hit_layToActiveRanks(lay, ranks)
Definition: hit_layout.h:1054
void hit_comStartSend(HitCom *issue)
Definition: hit_com.c:2051