Commits

John Wise committed fe040ef Merge

Merging

Comments (0)

Files changed (14)

src/enzo/CommunicationTransferParticles.C

-/***********************************************************************
-/
-/  COMMUNICATION ROUTINE: TRANSFER PARTICLES
-/
-/  written by: Greg Bryan
-/  date:       December, 1997
-/  modified:
-/
-/  PURPOSE:
-/
-************************************************************************/
-#ifndef OPTIMIZED_CTP 
-
-#ifdef USE_MPI
-#include "mpi.h"
-#endif /* USE_MPI */
- 
-#include <stdio.h>
-#include <string.h>
-#include <stdlib.h>
- 
-#include "ErrorExceptions.h"
-#include "macros_and_parameters.h"
-#include "typedefs.h"
-#include "global_data.h"
-#include "Fluxes.h"
-#include "GridList.h"
-#include "ExternalBoundary.h"
-#include "Grid.h"
-#include "TopGridData.h"
-#include "Hierarchy.h"
-#include "LevelHierarchy.h"
-void my_exit(int status);
- 
-// function prototypes
- 
-#ifdef USE_MPI
-static int FirstTimeCalled = TRUE;
-static MPI_Datatype MPI_ParticleMoveList;
-#endif
- 
- 
- 
- 
-int CommunicationTransferParticles(grid *GridPointer[], int NumberOfGrids)
-{
-
-  if (NumberOfGrids == 1)
-    return SUCCESS;
- 
-  struct ParticleMoveList {
-    int FromGrid;
-    int ToGrid[6];
-    int NumberToMove[6];
-    float_int *Pointer[6];
-  };
- 
-  /* This is a loop that keeps going until no particles were moved. */
- 
-  int NumberOfParticlesMoved, Done = FALSE;
- 
-  while (!Done) {
- 
-  NumberOfParticlesMoved = 0;
- 
-  int i, j, grid, GridsToSend = 0;
-  ParticleMoveList *SendList = new ParticleMoveList[NumberOfGrids];
- 
-  for (grid = 0; grid < NumberOfGrids; grid++)
-    for (i = 0; i < 6; i++) {
-      SendList[grid].NumberToMove[i] = 0;
-      SendList[grid].Pointer[i]      = NULL;
-      SendList[grid].ToGrid[i]       = -1;
-    }
-
-  /*
-  if ( MyProcessorNumber == 0 ) 
-  for (j = 0; j < NumberOfGrids; j++)
-    fprintf(stderr, "Pa(%"ISYM") CTP grid[%"ISYM"] = %"ISYM"\n",
-                    MyProcessorNumber, j, GridPointer[j]->ReturnNumberOfParticles());
-  */
- 
-  /* Generate the list of particle moves. */
- 
-  for (grid = 0; grid < NumberOfGrids; grid++)
-    if (GridPointer[grid]->ReturnProcessorNumber() == MyProcessorNumber) {
-      SendList[GridsToSend].FromGrid = grid;
-      if (GridPointer[grid]->CommunicationTransferParticles(GridPointer,
-	     NumberOfGrids,
-             SendList[GridsToSend].ToGrid,
-             SendList[GridsToSend].NumberToMove,
-             SendList[GridsToSend].Pointer,
-	     COPY_OUT) == FAIL) {
-	ENZO_FAIL("Error in grid->CommunicationTransferParticles.\n");
-      }
-      GridsToSend++;
-    }
- 
-  /* Share the Particle moves. */
- 
-  /* Allocate the array to receive subgrids. */
- 
-  ParticleMoveList *SharedList = NULL;
- 
-#ifdef USE_MPI
-  int NumberOfSharedGrids = 0;
-#endif /* USE_MPI */
-
-#ifdef USE_MPI
-  MPI_Status Status;
-  MPI_Datatype DataType = (sizeof(float) == 4) ? MPI_FLOAT : MPI_DOUBLE;
-  MPI_Datatype DataTypeInt = (sizeof(int) == 4) ? MPI_INT : MPI_LONG_LONG_INT;
-  MPI_Datatype DataTypeByte = MPI_BYTE;
-  MPI_Arg Count;
-  MPI_Arg SendCount;
-  MPI_Arg RecvCount;
-  MPI_Arg Source;
-  MPI_Arg Dest;
-  MPI_Arg Here;
-  MPI_Arg Usize;
-  MPI_Arg Xsize;
-  MPI_Arg Tag;
-  MPI_Arg stat;
-#endif /* USE_MPI */
- 
-  if (NumberOfProcessors > 1) {
- 
-#ifdef USE_MPI
- 
-    /* Generate a new MPI type corresponding to the ParticleMoveList struct. */
- 
-    if (FirstTimeCalled) {
-      Count = sizeof(ParticleMoveList);
-      //  fprintf(stderr, "Size of ParticleMoveList %"ISYM"\n", Count);
-      stat = MPI_Type_contiguous(Count, DataTypeByte, &MPI_ParticleMoveList);
-        if( stat != MPI_SUCCESS ){my_exit(EXIT_FAILURE);}
-      stat = MPI_Type_commit(&MPI_ParticleMoveList);
-        if( stat != MPI_SUCCESS ){my_exit(EXIT_FAILURE);}
-      FirstTimeCalled = FALSE;
-    }
-
-    int *SendListCount = new int[NumberOfProcessors];
-
-    MPI_Arg *MPI_SendListCount = new MPI_Arg[NumberOfProcessors];
-    MPI_Arg *MPI_SendListDisplacements = new MPI_Arg[NumberOfProcessors];
-
-    int jjj;
-
-    for ( jjj=0; jjj<NumberOfProcessors; jjj++) {
-      SendListCount[jjj] = 0;
-      MPI_SendListCount[jjj] = 0;
-      MPI_SendListDisplacements[jjj] = 0;
-    }
- 
-    /* Get counts from each processor to allocate buffers. */
- 
-#ifdef MPI_INSTRUMENTATION
-    starttime = MPI_Wtime();
-#endif /* MPI_INSTRUMENTATION */
-
-    SendCount = 1;
-    RecvCount = 1;
-
-    //fprintf(stderr, "GridsToSend %"ISYM" on %"ISYM"\n", GridsToSend, MyProcessorNumber);
-
-    stat = MPI_Allgather(&GridsToSend, SendCount, DataTypeInt, SendListCount, RecvCount, DataTypeInt, MPI_COMM_WORLD);
-      if( stat != MPI_SUCCESS ){my_exit(EXIT_FAILURE);}
-
-/*
-    if ( MyProcessorNumber == 0 )
-    for ( jjj=0; jjj<NumberOfProcessors; jjj++) {
-      fprintf(stderr, "SendListCount %"ISYM" on %"ISYM" \n", SendListCount[jjj], jjj);
-    } 
-*/
-
-    /* Allocate buffers and generated displacement list. */
- 
-    for (i = 0; i < NumberOfProcessors; i++) {
-      MPI_SendListDisplacements[i] = NumberOfSharedGrids;
-      NumberOfSharedGrids += SendListCount[i];
-      MPI_SendListCount[i] = SendListCount[i];
-    }
-    if (NumberOfSharedGrids != NumberOfGrids) {
-      ENZO_FAIL("CTP error\n");
-    }
-
-/*
-    if ( MyProcessorNumber == 0 )
-    for ( jjj=0; jjj<NumberOfProcessors; jjj++) {
-      fprintf(stderr, "SendListCount %"ISYM" SendListDisp %"ISYM" on %"ISYM" \n", MPI_SendListCount[jjj], MPI_SendListDisplacements[jjj], jjj);
-    }
-*/
-
-
-    SharedList = new ParticleMoveList[NumberOfSharedGrids];
- 
-    /* Perform sharing operation. */
-
-    Count = GridsToSend;
-     
-    stat = MPI_Allgatherv(SendList, Count, MPI_ParticleMoveList, SharedList,
-			  MPI_SendListCount, MPI_SendListDisplacements,
-			  MPI_ParticleMoveList, MPI_COMM_WORLD);
-      if( stat != MPI_SUCCESS ){my_exit(EXIT_FAILURE);}
- 
-#ifdef MPI_INSTRUMENTATION
-    endtime = MPI_Wtime();
-    timer[9] += endtime-starttime;
-    counter[9] ++;
-    timer[10] += double(NumberOfSharedGrids);
-    GlobalCommunication += endtime-starttime;
-    CommunicationTime += endtime-starttime;
-#endif /* MPI_INSTRUMENTATION */
- 
-    delete [] SendListCount;
-//    delete [] SendListDisplacements;
-    delete [] MPI_SendListCount;
-    delete [] MPI_SendListDisplacements;
- 
-#endif /* USE_MPI */
- 
-    /* Move particles. */
- 
-    for (j = 0; j < NumberOfGrids; j++) {
- 
-      int FromGrid = SharedList[j].FromGrid;
- 
-      /* Loop over transfer directions. */
- 
-      for (i = 0; i < 6; i++) {
- 
-	int ToGrid = max(SharedList[j].ToGrid[i], 0);
-	int FromProcessor = GridPointer[FromGrid]->ReturnProcessorNumber();
-	int ToProcessor = GridPointer[ToGrid]->ReturnProcessorNumber();
-	int TransferSize = SharedList[j].NumberToMove[i]*
-	                   (9+NumberOfParticleAttributes);
-
-	NumberOfParticlesMoved += SharedList[j].NumberToMove[i];
-
-	if (TransferSize > 0 && FromProcessor != ToProcessor) {
- 
-        //  fprintf(stderr, "P%"ISYM" I=%"ISYM" XFER SIZE %"ISYM"\n", MyProcessorNumber, i, TransferSize);
-
-#ifdef USE_MPI	
- 
-	  /* Send particles (Note: this is not good -- the data transfer will not work
-                across heterogeneous machines). */
- 
-#ifdef MPI_INSTRUMENTATION
-          starttime = MPI_Wtime();
-#endif  /* MPI_INSTRUMENTATION */
-
-	  Usize = sizeof(float_int);
-	  Xsize = TransferSize;
-
-          // fprintf(stderr, "sizeof SendCount %"ISYM"\n", sizeof(SendCount));
-          // fprintf(stderr, "sizeof MPI_Arg %"ISYM"\n", sizeof(MPI_Arg));
- 
-	  if (FromProcessor == MyProcessorNumber) {
-
-//          SendCount = TransferSize*sizeof(float_int);
-	    SendCount = Usize * Xsize;
-            Dest = ToProcessor;
-	    Here = MyProcessorNumber;
-
-	    if( SendCount > 0 ) {
-	    // fprintf(stderr, "P%"ISYM" Sending %"ISYM" bytes to %"ISYM" [%"ISYM" x %"ISYM"]\n", Here, SendCount, Dest, Usize, Xsize); 
-	    stat = MPI_Send(SharedList[j].Pointer[i], SendCount, DataTypeByte, Dest, MPI_TRANSFERPARTICLE_TAG, MPI_COMM_WORLD);
-	      if( stat != MPI_SUCCESS ){my_exit(EXIT_FAILURE);}
-	    }
-	  }
- 
-	  /* Receive particles. */
- 
-	  if (ToProcessor == MyProcessorNumber) {
-	    SharedList[j].Pointer[i] = new float_int[TransferSize];
-
-//          RecvCount = TransferSize*sizeof(float_int);
-	    RecvCount = Xsize*Usize;
-            Source = FromProcessor;
-	    Here = MyProcessorNumber;
-
-	    // fprintf(stderr, "P%"ISYM" Post receive %"ISYM" bytes from %"ISYM" [%"ISYM" x %"ISYM"]\n", Here, RecvCount, Source, Usize, Xsize);
-	    stat = MPI_Recv(SharedList[j].Pointer[i], RecvCount, DataTypeByte, Source, MPI_TRANSFERPARTICLE_TAG, MPI_COMM_WORLD, &Status);
-	      if( stat != MPI_SUCCESS ){my_exit(EXIT_FAILURE);}
-	    stat = MPI_Get_count(&Status, DataTypeByte, &RecvCount);
-	      if( stat != MPI_SUCCESS ){my_exit(EXIT_FAILURE);}
-	    // fprintf(stderr, "P%"ISYM" Received %"ISYM" bytes from %"ISYM" [%"ISYM"]\n", Here, RecvCount, Source, Usize);
-	  }
- 
-#ifdef MPI_INSTRUMENTATION
-          double endtime = MPI_Wtime();
-	  RecvComm += endtime-starttime;
-	  CommunicationTime += endtime-starttime;
-#endif  /* MPI_INSTRUMENTATION */
- 
-#endif /* USE_MPI */
- 
-	} // end: if (TransferSize > 0...
- 
-	/* If this is not on my processor, then clean up pointer, since it
-	   doesn't make sense on this processor (or if nothing moved). */
- 
-	if ((MyProcessorNumber != FromProcessor &&
-	     MyProcessorNumber != ToProcessor) || TransferSize == 0)
-	  SharedList[j].Pointer[i] = NULL;
- 
-      } // end: loop over i (directions)
-    } // end: loop over j (grids)
- 
-  } else {
-    SharedList = SendList;  // if there is only one processor
-    for (grid = 0; grid < NumberOfGrids; grid++)
-      for (i = 0; i < 6; i++)
-	NumberOfParticlesMoved += SharedList[grid].NumberToMove[i];
-  }
- 
-  /* Copy particles back to grids. */
- 
-  for (grid = 0; grid < NumberOfGrids; grid++)
-    if (GridPointer[grid]->ReturnProcessorNumber() == MyProcessorNumber) {
- 
-      int LocalNumberToMove[6], counter = 0;
-      float_int *LocalPointer[6];
-      for (i = 0; i < 6; i++)
-	LocalNumberToMove[i] = 0;
- 
-      /* Extract grid lists to put into this grid. */
- 
-      for (j = 0; j < NumberOfGrids; j++)
-	for (i = 0; i < 6; i++)
-	  if (SharedList[j].ToGrid[i] == grid) {
-	    LocalNumberToMove[counter] = SharedList[j].NumberToMove[i];
-	    LocalPointer[counter++] = SharedList[j].Pointer[i];
-	  }
- 
-      /* Copy particles back (SharedList[grid].ToGrid is a dummy, not used). */
-
-      if (GridPointer[grid]->CommunicationTransferParticles(GridPointer,
-	      NumberOfGrids, SharedList[grid].ToGrid, LocalNumberToMove,
-	      LocalPointer, COPY_IN) == FAIL) {
-	ENZO_FAIL("Error in grid->CommunicationTransferParticless\n");
-      }
- 
-    } // end: if grid is on my processor
- 
-  /* Set number of particles so everybody agrees. */
-
-#ifdef UNUSED 
-  if (NumberOfProcessors > 1) {
-    int *Changes = new int[NumberOfGrids];
-    for (j = 0; j < NumberOfGrids; j++)
-      Changes[j] = 0;
-    for (j = 0; j < NumberOfGrids; j++)
-      for (i = 0; i < 6; i++)
-	if (SharedList[j].ToGrid[i] != -1) {
-	  Changes[SharedList[j].FromGrid] -= SharedList[j].NumberToMove[i];
-	  Changes[SharedList[j].ToGrid[i]] += SharedList[j].NumberToMove[i];
-	}
-    for (j = 0; j < NumberOfGrids; j++) {
-      if (GridPointer[j]->ReturnProcessorNumber() != MyProcessorNumber)
-	GridPointer[j]->SetNumberOfParticles(
-		         GridPointer[j]->ReturnNumberOfParticles()+Changes[j]);
-      //      printf("Pb(%"ISYM") CTP grid[%"ISYM"] = %"ISYM"\n", MyProcessorNumber, j, GridPointer[j]->ReturnNumberOfParticles());
-    }
-    delete [] Changes;
-  }
-#endif
- 
-  /* CleanUp. */
- 
-  for (j = 0; j < NumberOfGrids; j++)
-    for (i = 0; i < 6; i++) {
-      delete [] SharedList[j].Pointer[i];
-    }
- 
-  if (SendList != SharedList)
-    delete [] SendList;
- 
-  delete [] SharedList;
- 
-  /* Check for completion. */
- 
-  if (debug)
-    printf("CommunicationTransferParticles: moved = %"ISYM"\n",
-	   NumberOfParticlesMoved);
-  if (NumberOfParticlesMoved == 0)
-
-    Done = TRUE;
- 
-  }
- 
-  return SUCCESS;
-}
-#endif /* OPTIMIZED_CTP */

src/enzo/CommunicationTransferParticlesOpt.C

 /  PURPOSE:
 /
 ************************************************************************/
-#ifdef OPTIMIZED_CTP
- 
 #ifdef USE_MPI
 #include "mpi.h"
 #endif /* USE_MPI */
 #include <stdio.h>
 #include <string.h>
 #include <stdlib.h>
+//#include <algorithm>
  
 #include "ErrorExceptions.h"
 #include "macros_and_parameters.h"
 #include "Hierarchy.h"
 #include "LevelHierarchy.h"
 #include "CommunicationUtilities.h"
+
 void my_exit(int status);
  
 // function prototypes
 int CommunicationShareParticles(int *NumberToMove, particle_data* &SendList,
 				int &NumberOfReceives,
 				particle_data* &SharedList);
+int search_lower_bound(int *arr, int value, int low, int high, 
+		       int total);
 
 #define NO_DEBUG_CTP
 #define KEEP_PARTICLES_LOCAL
  
-int CommunicationTransferParticles(grid *GridPointer[], int NumberOfGrids)
+int CommunicationTransferParticles(grid *GridPointer[], int NumberOfGrids,
+				   int TopGridDims[])
 {
 
   if (NumberOfGrids == 1)
     return SUCCESS;
  
-  int i, j, jstart, jend, dim, grid, proc;
+  int i, j, jstart, jend, dim, grid, proc, DisplacementCount, ThisCount;
+  float ExactDims, ExactCount;
 
   /* Assume that the grid was split by Enzo_Dims_create, and create a
      map from grid number to an index that is determined from (i,j,k)
      of the grid partitions. */
 
   int *GridMap = new int[NumberOfGrids];
+  int *StartIndex[MAX_DIMENSION];
   int Layout[MAX_DIMENSION], LayoutTemp[MAX_DIMENSION];
-  int Rank, grid_num, GridPosition[MAX_DIMENSION], Dims[MAX_DIMENSION];
+  int Rank, grid_num, CenterIndex, bin;
+  int *pbin;
+  int GridPosition[MAX_DIMENSION], Dims[MAX_DIMENSION];
   FLOAT Left[MAX_DIMENSION], Right[MAX_DIMENSION];
 
   for (dim = 0; dim < MAX_DIMENSION; dim++) {
   for (dim = Rank; dim < MAX_DIMENSION; dim++)
     Layout[Rank-1-dim] = 0;
 
+  /* For unequal splits of the topgrid, we need the start indices of
+     the partitions.  It's easier to recalculate than to search for
+     them. */
+
+  for (dim = 0; dim < Rank; dim++) {
+
+    StartIndex[dim] = new int[Layout[dim]+1];
+    ExactDims = float(TopGridDims[dim]) / float(Layout[dim]);
+    ExactCount = 0.0;
+    DisplacementCount = 0;
+
+    for (i = 0; i < Layout[dim]; i++) {
+      ExactCount += ExactDims;
+      if (dim == 0)
+	ThisCount = nint(0.5*ExactCount)*2 - DisplacementCount;
+      else
+	ThisCount = nint(ExactCount) - DisplacementCount;
+      StartIndex[dim][i] = DisplacementCount;
+      DisplacementCount += ThisCount;
+    } // ENDFOR i
+
+    StartIndex[dim][Layout[dim]] = TopGridDims[dim];
+
+  } // ENDFOR dim
+
   for (grid = 0; grid < NumberOfGrids; grid++) {
     GridPointer[grid]->ReturnGridInfo(&Rank, Dims, Left, Right);
-    for (dim = 0; dim < Rank; dim++)
-      GridPosition[dim] = 
-	int(Layout[dim] * (0.5*(Right[dim]+Left[dim]) - DomainLeftEdge[dim]) /
-	    (DomainRightEdge[dim] - DomainLeftEdge[dim]));
+    for (dim = 0; dim < Rank; dim++) {
+
+      if (Layout[dim] == 1) {
+	GridPosition[dim] = 0;
+      } else {
+
+	CenterIndex = 
+	  int(TopGridDims[dim] *
+	      (0.5*(Right[dim]+Left[dim]) - DomainLeftEdge[dim]) /
+	      (DomainRightEdge[dim] - DomainLeftEdge[dim]));
+      
+	GridPosition[dim] = 
+	  search_lower_bound(StartIndex[dim], CenterIndex, 0, Layout[dim],
+			     Layout[dim]);
+
+      } // ENDELSE
+
+    } // ENDFOR dim
+//    if (debug)
+//      printf("grid %d: GPos = %d %d %d\n", grid, GridPosition[0],
+//	     GridPosition[1], GridPosition[2]);
     grid_num = GridPosition[0] + 
       Layout[0] * (GridPosition[1] + Layout[1]*GridPosition[2]);
     GridMap[grid_num] = grid;
   for (grid = 0; grid < NumberOfGrids; grid++)
     if (GridPointer[grid]->
 	CommunicationTransferParticles(GridPointer, NumberOfGrids, grid, 
-				       NumberToMove, Zero, Zero, SendList, 
-				       Layout, GridMap, COPY_OUT) == FAIL) {
+				       TopGridDims, NumberToMove, Zero, 
+				       Zero, SendList, Layout, StartIndex, 
+				       GridMap, COPY_OUT) == FAIL) {
       ENZO_FAIL("Error in grid->CommunicationTransferParticles(OUT).\n");
     }
 
   NumberOfReceives = TotalNumberToMove;
   int particle_data_size = sizeof(particle_data);
   qsort(SharedList, TotalNumberToMove, particle_data_size, compare_grid);
+  //std::sort(SharedList, SharedList+TotalNumberToMove, cmp_grid());
 
 #else
 
 	if (jend == NumberOfReceives) break;
       }
       if (GridPointer[j]->CommunicationTransferParticles
-	  (GridPointer, NumberOfGrids, j, NumberToMove, jstart, jend, 
-	   SharedList, Layout, GridMap, COPY_IN) == FAIL) {
+	  (GridPointer, NumberOfGrids, j, TopGridDims,
+	   NumberToMove, jstart, jend, 
+	   SharedList, Layout, StartIndex, GridMap, COPY_IN) == FAIL) {
 	ENZO_FAIL("Error in grid->CommunicationTransferParticles(IN).\n");
       }
       jstart = jend;
   delete [] SharedList;
   delete [] NumberToMove;
   delete [] GridMap;
+  for (dim = 0; dim < Rank; dim++)
+    delete [] StartIndex[dim];
 
   CommunicationSumValues(&TotalNumberToMove, 1);
   if (debug)
   return SUCCESS;
 }
 
-#endif /* OPTIMIZED_CTP */

src/enzo/CommunicationTransferStars.C

-/***********************************************************************
-/
-/  COMMUNICATION ROUTINE: TRANSFER STARS
-/
-/  written by: Greg Bryan
-/  date:       December, 1997
-/  modified:   John Wise
-/  date:       July, 2009 (modified CTP)
-/
-/  PURPOSE:
-/
-************************************************************************/
-#ifndef OPTIMIZED_CTP 
-
-#ifdef USE_MPI
-#include "mpi.h"
-#endif /* USE_MPI */
- 
-#include <stdio.h>
-#include <string.h>
-#include <stdlib.h>
- 
-#include "ErrorExceptions.h"
-#include "macros_and_parameters.h"
-#include "typedefs.h"
-#include "global_data.h"
-#include "Fluxes.h"
-#include "GridList.h"
-#include "ExternalBoundary.h"
-#include "Grid.h"
-#include "TopGridData.h"
-#include "Hierarchy.h"
-#include "LevelHierarchy.h"
-void my_exit(int status);
- 
-// function prototypes
- 
-#ifdef USE_MPI
-static int FirstTimeCalled = TRUE;
-static MPI_Datatype MPI_StarMoveList;
-#endif
- 
- 
- 
- 
-int CommunicationTransferStars(grid *GridPointer[], int NumberOfGrids)
-{
-
-  if (NumberOfGrids == 1)
-    return SUCCESS;
- 
-  struct StarMoveList {
-    int FromGrid;
-    int ToGrid[6];
-    int NumberToMove[6];
-    StarBuffer *Pointer[6];
-  };
- 
-  /* This is a loop that keeps going until no stars were moved. */
- 
-  int NumberOfStarsMoved, Done = FALSE;
- 
-  while (!Done) {
- 
-  NumberOfStarsMoved = 0;
- 
-  int i, j, grid, GridsToSend = 0;
-  StarMoveList *SendList = new StarMoveList[NumberOfGrids];
- 
-  for (grid = 0; grid < NumberOfGrids; grid++)
-    for (i = 0; i < 6; i++) {
-      SendList[grid].NumberToMove[i] = 0;
-      SendList[grid].Pointer[i]      = NULL;
-      SendList[grid].ToGrid[i]       = -1;
-    }
-
-  /*
-  if ( MyProcessorNumber == 0 ) 
-  for (j = 0; j < NumberOfGrids; j++)
-    fprintf(stderr, "Pa(%"ISYM") CTP grid[%"ISYM"] = %"ISYM"\n",
-                    MyProcessorNumber, j, GridPointer[j]->ReturnNumberOfStars());
-  */
- 
-  /* Generate the list of star moves. */
- 
-  for (grid = 0; grid < NumberOfGrids; grid++)
-    if (GridPointer[grid]->ReturnProcessorNumber() == MyProcessorNumber) {
-      SendList[GridsToSend].FromGrid = grid;
-      GridPointer[grid]->CommunicationTransferStars
-	(GridPointer, NumberOfGrids, 
-	 SendList[GridsToSend].ToGrid, 
-	 SendList[GridsToSend].NumberToMove, 
-	 SendList[GridsToSend].Pointer, COPY_OUT);
-      GridsToSend++;
-    }
- 
-  /* Share the star moves. */
- 
-  /* Allocate the array to receive subgrids. */
- 
-  StarMoveList *SharedList = NULL;
- 
-#ifdef USE_MPI
-  int NumberOfSharedGrids = 0;
-#endif /* USE_MPI */
-
-#ifdef USE_MPI
-  MPI_Status Status;
-  MPI_Datatype DataType = (sizeof(float) == 4) ? MPI_FLOAT : MPI_DOUBLE;
-  MPI_Datatype DataTypeInt = (sizeof(int) == 4) ? MPI_INT : MPI_LONG_LONG_INT;
-  MPI_Datatype DataTypeByte = MPI_BYTE;
-  MPI_Arg Count;
-  MPI_Arg SendCount;
-  MPI_Arg RecvCount;
-  MPI_Arg Source;
-  MPI_Arg Dest;
-  MPI_Arg Here;
-  MPI_Arg Usize;
-  MPI_Arg Xsize;
-  MPI_Arg Tag;
-  MPI_Arg stat;
-#endif /* USE_MPI */
- 
-  if (NumberOfProcessors > 1) {
- 
-#ifdef USE_MPI
- 
-    /* Generate a new MPI type corresponding to the StarMoveList struct. */
- 
-    if (FirstTimeCalled) {
-      Count = sizeof(StarMoveList);
-      //  fprintf(stderr, "Size of StarMoveList %"ISYM"\n", Count);
-      stat = MPI_Type_contiguous(Count, DataTypeByte, &MPI_StarMoveList);
-        if( stat != MPI_SUCCESS ){my_exit(EXIT_FAILURE);}
-      stat = MPI_Type_commit(&MPI_StarMoveList);
-        if( stat != MPI_SUCCESS ){my_exit(EXIT_FAILURE);}
-      FirstTimeCalled = FALSE;
-    }
-
-    int *SendListCount = new int[NumberOfProcessors];
-
-    MPI_Arg *MPI_SendListCount = new MPI_Arg[NumberOfProcessors];
-    MPI_Arg *MPI_SendListDisplacements = new MPI_Arg[NumberOfProcessors];
-
-    int jjj;
-
-    for ( jjj=0; jjj<NumberOfProcessors; jjj++) {
-      SendListCount[jjj] = 0;
-      MPI_SendListCount[jjj] = 0;
-      MPI_SendListDisplacements[jjj] = 0;
-    }
- 
-    /* Get counts from each processor to allocate buffers. */
- 
-#ifdef MPI_INSTRUMENTATION
-    starttime = MPI_Wtime();
-#endif /* MPI_INSTRUMENTATION */
-
-    SendCount = 1;
-    RecvCount = 1;
-
-    //fprintf(stderr, "GridsToSend %"ISYM" on %"ISYM"\n", GridsToSend, MyProcessorNumber);
-
-    stat = MPI_Allgather(&GridsToSend, SendCount, DataTypeInt, SendListCount, RecvCount, DataTypeInt, MPI_COMM_WORLD);
-      if( stat != MPI_SUCCESS ){my_exit(EXIT_FAILURE);}
-
-/*
-    if ( MyProcessorNumber == 0 )
-    for ( jjj=0; jjj<NumberOfProcessors; jjj++) {
-      fprintf(stderr, "SendListCount %"ISYM" on %"ISYM" \n", SendListCount[jjj], jjj);
-    } 
-*/
-
-    /* Allocate buffers and generated displacement list. */
- 
-    for (i = 0; i < NumberOfProcessors; i++) {
-      MPI_SendListDisplacements[i] = NumberOfSharedGrids;
-      NumberOfSharedGrids += SendListCount[i];
-      MPI_SendListCount[i] = SendListCount[i];
-    }
-    if (NumberOfSharedGrids != NumberOfGrids) {
-      ENZO_FAIL("CTP error\n");
-    }
-
-/*
-    if ( MyProcessorNumber == 0 )
-    for ( jjj=0; jjj<NumberOfProcessors; jjj++) {
-      fprintf(stderr, "SendListCount %"ISYM" SendListDisp %"ISYM" on %"ISYM" \n", MPI_SendListCount[jjj], MPI_SendListDisplacements[jjj], jjj);
-    }
-*/
-
-
-    SharedList = new StarMoveList[NumberOfSharedGrids];
- 
-    /* Perform sharing operation. */
-
-    Count = GridsToSend;
-     
-    stat = MPI_Allgatherv(SendList, Count, MPI_StarMoveList, SharedList,
-			  MPI_SendListCount, MPI_SendListDisplacements,
-			  MPI_StarMoveList, MPI_COMM_WORLD);
-      if( stat != MPI_SUCCESS ){my_exit(EXIT_FAILURE);}
- 
-#ifdef MPI_INSTRUMENTATION
-    endtime = MPI_Wtime();
-    timer[9] += endtime-starttime;
-    counter[9] ++;
-    timer[10] += double(NumberOfSharedGrids);
-    GlobalCommunication += endtime-starttime;
-    CommunicationTime += endtime-starttime;
-#endif /* MPI_INSTRUMENTATION */
- 
-    delete [] SendListCount;
-//    delete [] SendListDisplacements;
-    delete [] MPI_SendListCount;
-    delete [] MPI_SendListDisplacements;
- 
-#endif /* USE_MPI */
- 
-    /* Move stars. */
- 
-    for (j = 0; j < NumberOfGrids; j++) {
- 
-      int FromGrid = SharedList[j].FromGrid;
- 
-      /* Loop over transfer directions. */
- 
-      for (i = 0; i < 6; i++) {
- 
-	int ToGrid = max(SharedList[j].ToGrid[i], 0);
-	int FromProcessor = GridPointer[FromGrid]->ReturnProcessorNumber();
-	int ToProcessor = GridPointer[ToGrid]->ReturnProcessorNumber();
-	int TransferSize = SharedList[j].NumberToMove[i];
-
-	NumberOfStarsMoved += SharedList[j].NumberToMove[i];
-
-	if (TransferSize > 0 && FromProcessor != ToProcessor) {
- 
-        //  fprintf(stderr, "P%"ISYM" I=%"ISYM" XFER SIZE %"ISYM"\n", MyProcessorNumber, i, TransferSize);
-
-#ifdef USE_MPI	
- 
-	  /* Send stars (Note: this is not good -- the data transfer
-	     will not work across heterogeneous machines). */
- 
-#ifdef MPI_INSTRUMENTATION
-          starttime = MPI_Wtime();
-#endif  /* MPI_INSTRUMENTATION */
-
-	  Usize = sizeof(StarBuffer);
-	  Xsize = TransferSize;
-
-          // fprintf(stderr, "sizeof SendCount %"ISYM"\n", sizeof(SendCount));
-          // fprintf(stderr, "sizeof MPI_Arg %"ISYM"\n", sizeof(MPI_Arg));
- 
-	  if (FromProcessor == MyProcessorNumber) {
-
-//          SendCount = TransferSize*sizeof(StarBuffer);
-	    SendCount = Usize * Xsize;
-            Dest = ToProcessor;
-	    Here = MyProcessorNumber;
-
-	    if( SendCount > 0 ) {
-	    // fprintf(stderr, "P%"ISYM" Sending %"ISYM" bytes to %"ISYM" [%"ISYM" x %"ISYM"]\n", Here, SendCount, Dest, Usize, Xsize); 
-	    stat = MPI_Send(SharedList[j].Pointer[i], SendCount, DataTypeByte, Dest, MPI_TRANSFERPARTICLE_TAG, MPI_COMM_WORLD);
-	      if( stat != MPI_SUCCESS ){my_exit(EXIT_FAILURE);}
-	    }
-	  }
- 
-	  /* Receive stars. */
- 
-	  if (ToProcessor == MyProcessorNumber) {
-	    SharedList[j].Pointer[i] = new StarBuffer[TransferSize];
-
-//          RecvCount = TransferSize*sizeof(StarBuffer);
-	    RecvCount = Xsize*Usize;
-            Source = FromProcessor;
-	    Here = MyProcessorNumber;
-
-	    // fprintf(stderr, "P%"ISYM" Post receive %"ISYM" bytes from %"ISYM" [%"ISYM" x %"ISYM"]\n", Here, RecvCount, Source, Usize, Xsize);
-	    stat = MPI_Recv(SharedList[j].Pointer[i], RecvCount, DataTypeByte, Source, MPI_TRANSFERPARTICLE_TAG, MPI_COMM_WORLD, &Status);
-	      if( stat != MPI_SUCCESS ){my_exit(EXIT_FAILURE);}
-	    stat = MPI_Get_count(&Status, DataTypeByte, &RecvCount);
-	      if( stat != MPI_SUCCESS ){my_exit(EXIT_FAILURE);}
-	    // fprintf(stderr, "P%"ISYM" Received %"ISYM" bytes from %"ISYM" [%"ISYM"]\n", Here, RecvCount, Source, Usize);
-	  }
- 
-#ifdef MPI_INSTRUMENTATION
-          double endtime = MPI_Wtime();
-	  RecvComm += endtime-starttime;
-	  CommunicationTime += endtime-starttime;
-#endif  /* MPI_INSTRUMENTATION */
- 
-#endif /* USE_MPI */
- 
-	} // end: if (TransferSize > 0...
- 
-	/* If this is not on my processor, then clean up pointer, since it
-	   doesn't make sense on this processor (or if nothing moved). */
- 
-	if ((MyProcessorNumber != FromProcessor &&
-	     MyProcessorNumber != ToProcessor) || TransferSize == 0)
-	  SharedList[j].Pointer[i] = NULL;
- 
-      } // end: loop over i (directions)
-    } // end: loop over j (grids)
- 
-  } else {
-    SharedList = SendList;  // if there is only one processor
-    for (grid = 0; grid < NumberOfGrids; grid++)
-      for (i = 0; i < 6; i++)
-	NumberOfStarsMoved += SharedList[grid].NumberToMove[i];
-  }
- 
-  /* Copy stars back to grids. */
- 
-  for (grid = 0; grid < NumberOfGrids; grid++)
-    if (GridPointer[grid]->ReturnProcessorNumber() == MyProcessorNumber) {
- 
-      int LocalNumberToMove[6], counter = 0;
-      StarBuffer *LocalPointer[6];
-      for (i = 0; i < 6; i++)
-	LocalNumberToMove[i] = 0;
- 
-      /* Extract grid lists to put into this grid. */
- 
-      for (j = 0; j < NumberOfGrids; j++)
-	for (i = 0; i < 6; i++)
-	  if (SharedList[j].ToGrid[i] == grid) {
-	    LocalNumberToMove[counter] = SharedList[j].NumberToMove[i];
-	    LocalPointer[counter++] = SharedList[j].Pointer[i];
-	  }
- 
-      /* Copy stars back (SharedList[grid].ToGrid is a dummy, not used). */
- 
-      GridPointer[grid]->CommunicationTransferStars
-	(GridPointer, NumberOfGrids, SharedList[grid].ToGrid, 
-	 LocalNumberToMove, LocalPointer, COPY_IN);
- 
-    } // end: if grid is on my processor
-
-  /* Set number of particles so everybody agrees. */
-
-#ifdef UNUSED 
-  if (NumberOfProcessors > 1) {
-    int *Changes = new int[NumberOfGrids];
-    for (j = 0; j < NumberOfGrids; j++)
-      Changes[j] = 0;
-    for (j = 0; j < NumberOfGrids; j++)
-      for (i = 0; i < 6; i++)
-	if (SharedList[j].ToGrid[i] != -1) {
-	  Changes[SharedList[j].FromGrid] -= SharedList[j].NumberToMove[i];
-	  Changes[SharedList[j].ToGrid[i]] += SharedList[j].NumberToMove[i];
-	}
-    for (j = 0; j < NumberOfGrids; j++) {
-      if (GridPointer[j]->ReturnProcessorNumber() != MyProcessorNumber)
-	GridPointer[j]->SetNumberOfStars(
-		         GridPointer[j]->ReturnNumberOfStars()+Changes[j]);
-      //      printf("Pb(%"ISYM") CTP grid[%"ISYM"] = %"ISYM"\n", MyProcessorNumber, j, GridPointer[j]->ReturnNumberOfStars());
-    }
-    delete [] Changes;
-  }
-#endif /* UNUSED */
- 
-  /* Set number of stars so everybody agrees. */
- 
-  if (NumberOfProcessors > 1) {
-    int *Changes = new int[NumberOfGrids];
-    for (j = 0; j < NumberOfGrids; j++)
-      Changes[j] = 0;
-    for (j = 0; j < NumberOfGrids; j++)
-      for (i = 0; i < 6; i++)
-	if (SharedList[j].ToGrid[i] != -1) {
-	  Changes[SharedList[j].FromGrid] -= SharedList[j].NumberToMove[i];
-	  Changes[SharedList[j].ToGrid[i]] += SharedList[j].NumberToMove[i];
-	}
-    for (j = 0; j < NumberOfGrids; j++) {
-      if (GridPointer[j]->ReturnProcessorNumber() != MyProcessorNumber)
-	GridPointer[j]->SetNumberOfStars(
-		         GridPointer[j]->ReturnNumberOfStars()+Changes[j]);
-      //      printf("Pb(%"ISYM") CTP grid[%"ISYM"] = %"ISYM"\n", MyProcessorNumber, j, GridPointer[j]->ReturnNumberOfParticles());
-    }
-    delete [] Changes;
-  }
-
-  /* CleanUp. */
- 
-  for (j = 0; j < NumberOfGrids; j++)
-    for (i = 0; i < 6; i++) {
-      delete [] SharedList[j].Pointer[i];
-    }
- 
-  if (SendList != SharedList)
-    delete [] SendList;
- 
-  delete [] SharedList;
- 
-  /* Check for completion. */
- 
-  if (debug)
-    printf("CommunicationTransferStars: moved = %"ISYM"\n",
-	   NumberOfStarsMoved);
-  if (NumberOfStarsMoved == 0)
-
-    Done = TRUE;
- 
-  }
- 
-  return SUCCESS;
-}
-#endif /* OPTIMIZED_CTP */

src/enzo/CommunicationTransferStarsOpt.C

 /  PURPOSE:
 /
 ************************************************************************/
-#ifdef OPTIMIZED_CTP
- 
 #ifdef USE_MPI
 #include "mpi.h"
 #endif /* USE_MPI */
 #include <stdio.h>
 #include <string.h>
 #include <stdlib.h>
+//#include <algorithm>
  
 #include "ErrorExceptions.h"
 #include "macros_and_parameters.h"
 #include "Hierarchy.h"
 #include "LevelHierarchy.h"
 #include "CommunicationUtilities.h"
+
 void my_exit(int status);
  
 // function prototypes
  
 Eint32 compare_star_grid(const void *a, const void *b);
 int Enzo_Dims_create(int nnodes, int ndims, int *dims);
+int search_lower_bound(int *arr, int value, int low, int high, 
+		       int total);
 
 // Remove define.  This method will always be used.
 //#define KEEP_PARTICLES_LOCAL
  
-int CommunicationTransferStars(grid *GridPointer[], int NumberOfGrids)
+int CommunicationTransferStars(grid *GridPointer[], int NumberOfGrids,
+			       int TopGridDims[])
 {
 
   if (NumberOfGrids == 1)
     return SUCCESS;
- 
-  int i, j, jstart, jend, dim, grid, proc;
+
+  int i, j, jstart, jend, dim, grid, proc, DisplacementCount, ThisCount;
+  float ExactDims, ExactCount;
 
   /* Assume that the grid was split by Enzo_Dims_create, and create a
      map from grid number to an index that is determined from (i,j,k)
      of the grid partitions. */
 
   int *GridMap = new int[NumberOfGrids];
+  int *StartIndex[MAX_DIMENSION];
   int Layout[MAX_DIMENSION], LayoutTemp[MAX_DIMENSION];
-  int Rank, grid_num, GridPosition[MAX_DIMENSION], Dims[MAX_DIMENSION];
+  int Rank, grid_num, bin, CenterIndex;
+  int *pbin;
+  int GridPosition[MAX_DIMENSION], Dims[MAX_DIMENSION];
   FLOAT Left[MAX_DIMENSION], Right[MAX_DIMENSION];
 
   for (dim = 0; dim < MAX_DIMENSION; dim++) {
     Layout[dim] = 0;
     GridPosition[dim] = 0;
+    StartIndex[dim] = NULL;
   }
 
   GridPointer[0]->ReturnGridInfo(&Rank, Dims, Left, Right); // need rank
   for (dim = Rank; dim < MAX_DIMENSION; dim++)
     Layout[Rank-1-dim] = 0;
 
+  /* For unequal splits of the topgrid, we need the start indices of
+     the partitions.  It's easier to recalculate than to search for
+     them. */
+
+  for (dim = 0; dim < Rank; dim++) {
+
+    StartIndex[dim] = new int[Layout[dim]+1];
+    ExactDims = float(TopGridDims[dim]) / float(Layout[dim]);
+    ExactCount = 0.0;
+    DisplacementCount = 0;
+
+    for (i = 0; i < Layout[dim]; i++) {
+      ExactCount += ExactDims;
+      if (dim == 0)
+	ThisCount = nint(0.5*ExactCount)*2 - DisplacementCount;
+      else
+	ThisCount = nint(ExactCount) - DisplacementCount;
+      StartIndex[dim][i] = DisplacementCount;
+      DisplacementCount += ThisCount;
+    } // ENDFOR i    
+
+    StartIndex[dim][Layout[dim]] = TopGridDims[dim];
+
+  } // ENDFOR dim
+
   for (grid = 0; grid < NumberOfGrids; grid++) {
     GridPointer[grid]->ReturnGridInfo(&Rank, Dims, Left, Right);
-    for (dim = 0; dim < Rank; dim++)
-      GridPosition[dim] = 
-	int(Layout[dim] * (0.5*(Right[dim]+Left[dim]) - DomainLeftEdge[dim]) /
-	    (DomainRightEdge[dim] - DomainLeftEdge[dim]));
+    for (dim = 0; dim < Rank; dim++) {
+
+      if (Layout[dim] == 1) {
+	GridPosition[dim] = 0;
+      } else {
+
+	CenterIndex = 
+	  int(TopGridDims[dim] *
+	      (0.5*(Right[dim]+Left[dim]) - DomainLeftEdge[dim]) /
+	      (DomainRightEdge[dim] - DomainLeftEdge[dim]));
+
+	GridPosition[dim] = 
+	  search_lower_bound(StartIndex[dim], CenterIndex, 0, Layout[dim],
+			     Layout[dim]);
+
+      } // ENDELSE
+
+    } // ENDFOR dim
     grid_num = GridPosition[0] + 
       Layout[0] * (GridPosition[1] + Layout[1]*GridPosition[2]);
     GridMap[grid_num] = grid;
   int Zero = 0;
   for (grid = 0; grid < NumberOfGrids; grid++)
     GridPointer[grid]->CommunicationTransferStars
-      (GridPointer, NumberOfGrids, grid, NumberToMove, Zero, Zero, 
-       SendList, Layout, GridMap, COPY_OUT);
+      (GridPointer, NumberOfGrids, grid, TopGridDims, NumberToMove, Zero, Zero, 
+       SendList, Layout, StartIndex, GridMap, COPY_OUT);
 
   int TotalNumberToMove = 0;
   for (proc = 0; proc < NumberOfProcessors; proc++)
   NumberOfReceives = TotalNumberToMove;
   int star_data_size = sizeof(star_data);
   qsort(SharedList, TotalNumberToMove, star_data_size, compare_star_grid);
+  //std::sort(SharedList, SharedList+TotalNumberToMove, cmp_star_grid());
 
   /* Copy stars back to grids */
 
       }
 
       GridPointer[j]->CommunicationTransferStars
-	(GridPointer, NumberOfGrids, j, NumberToMove, jstart, jend, 
-	 SharedList, Layout, GridMap, COPY_IN);
+	(GridPointer, NumberOfGrids, j, TopGridDims, NumberToMove, 
+	 jstart, jend, SharedList, Layout, StartIndex, GridMap, COPY_IN);
 
       jstart = jend;
     } // ENDFOR grids
   delete [] SharedList;
   delete [] NumberToMove;
   delete [] GridMap;
+  for (dim = 0; dim < Rank; dim++)
+    delete [] StartIndex[dim];
 
   CommunicationSumValues(&TotalNumberToMove, 1);
   if (debug)
     printf("CommunicationTransferStars: moved = %"ISYM"\n",
   	   TotalNumberToMove);
- 
+
   return SUCCESS;
 }
 
-#endif /* OPTIMIZED_CTP */

src/enzo/EvolveHierarchy.C

   MetaData.StartCPUTime = MetaData.CPUTime = LastCPUTime = ReturnWallTime();
   MetaData.LastCycleCPUTime = 0.0;
  
-  /* Double-check if the topgrid is evenly divided if we're using the
-     optimized version of CommunicationTransferParticles. */
-
-#ifdef OPTIMIZED_CTP
-  int NumberOfGrids = 0, Layout[MAX_DIMENSION] = {1,1,1};
-  Temp = LevelArray[0];
-  while (Temp != NULL) {
-    NumberOfGrids++;
-    Temp = Temp->NextGridThisLevel;
-  }
-  Enzo_Dims_create(NumberOfGrids, MetaData.TopGridRank, Layout);
-  for (dim = 0; dim < MetaData.TopGridRank; dim++)
-    if (MetaData.TopGridDims[dim] % Layout[MAX_DIMENSION-1-dim] != 0) {
-      if (debug)
-	fprintf(stderr, "ERROR: "
-		"\tTo use the optimized CommunicationTransferParticles routine,\n"
-		"\tthe top grid must be split evenly, "
-		"i.e. mod(Dims[i], Layout[i]) must be 0\n"
-		"\t==> dimension %"ISYM": Dims = %"ISYM", Layout = %"ISYM"\n",
-		dim, MetaData.TopGridDims[dim], Layout[MAX_DIMENSION-1-dim]);
-      ENZO_FAIL("");
-    }
-#endif /* OPTIMIZED_CTP */
-
   /* Attach RandomForcingFields to BaryonFields temporarily to apply BCs */
  
   if (RandomForcing) { //AK

src/enzo/FOF_Finalize.C

 int GetUnits(float *DensityUnits, float *LengthUnits,
 	     float *TemperatureUnits, float *TimeUnits,
 	     float *VelocityUnits, FLOAT Time);
-int CommunicationTransferParticles(grid *GridPointer[], int NumberOfGrids);
+int CommunicationTransferParticles(grid *GridPointer[], int NumberOfGrids,
+				   int TopGridDims[]);
 int CommunicationCollectParticles(LevelHierarchyEntry *LevelArray[],
 				  int level, bool ParticlesAreLocal, 
 				  bool SyncNumberOfParticles, 
 				SyncNumberOfParticles, MoveStars,
 				SIBLINGS_ONLY);
 
-  CommunicationTransferParticles(GridPointer, ngrids);
+  CommunicationTransferParticles(GridPointer, ngrids, MetaData->TopGridDims);
 
   ParticlesAreLocal = false;
   SyncNumberOfParticles = true;
 
 /* Transfer particle amount level 0 grids. */
 
-#ifdef OPTIMIZED_CTP
   int CommunicationTransferParticles(grid* Grids[], int NumberOfGrids,
-				     int ThisGridNum, int *&NumberToMove, 
+				     int ThisGridNum, int TopGridDims[],
+				     int *&NumberToMove, 
 				     int StartIndex, int EndIndex, 
 				     particle_data *&List,
-				     int *Layout, int *GridMap, 
-				     int CopyDirection);
+				     int *Layout, int *GStartIndex[],
+				     int *GridMap, int CopyDirection);
   int CommunicationTransferStars(grid* Grids[], int NumberOfGrids,
-				 int ThisGridNum, int *&NumberToMove, 
+				 int ThisGridNum, int TopGridDims[],
+				 int *&NumberToMove, 
 				 int StartIndex, int EndIndex, 
-				 star_data *&List,
-				 int *Layout, int *GridMap, 
+				 star_data *&List, int *Layout, 
+				 int *GStartIndex[], int *GridMap, 
 				 int CopyDirection);
-#else
-  int CommunicationTransferParticles(grid* Grids[], int NumberOfGrids,
-			     int ToGrid[6], int NumberToMove[6],
-			     float_int *ParticleData[6], int CopyDirection);
-  int CommunicationTransferStars(grid* Grids[], int NumberOfGrids,
-			 int ToGrid[6], int NumberToMove[6],
-			 StarBuffer *StarData[6], int CopyDirection);
-#endif
 
   int CollectParticles(int GridNum, int* &NumberToMove, 
 		       int &StartIndex, int &EndIndex, 

src/enzo/Grid_CommunicationTransferParticles.C

-/***********************************************************************
-/
-/  GRID CLASS (COPY PARTICLES INTO OR OUT OF GRID)
-/
-/  written by: Greg Bryan
-/  date:       January, 1999
-/  modified1:  Robert Harkness
-/  date:       April, 2006
-/
-/  PURPOSE:
-/
-************************************************************************/
-#ifndef OPTIMIZED_CTP 
-//
- 
-#ifdef USE_MPI
-#include <mpi.h>
-#endif
-#include <stdio.h>
-
-#include "ErrorExceptions.h"
-#include "macros_and_parameters.h"
-#include "typedefs.h"
-#include "global_data.h"
-#include "Fluxes.h"
-#include "GridList.h"
-#include "ExternalBoundary.h"
-#include "Grid.h"
-#include "Hierarchy.h"
- 
- 
-int grid::CommunicationTransferParticles(grid* Grids[], int NumberOfGrids,
-		 int ToGrid[6], int NumberToMove[6],
-		 float_int *ParticleData[6], int CopyDirection)
-{
- 
-  /* Declarations. */
- 
-  int i, j, k, dim, grid;
- 
-  /* ----------------------------------------------------------------- */
-  /* Copy particle out of grid. */
- 
-  if (CopyDirection == COPY_OUT) {
-
-    /* If there are no particles to move, we're done. */
- 
-    if (NumberOfParticles == 0)
-      return SUCCESS;
- 
-    /* Count particles to move (mark already counted by setting mass < 0). */
- 
-    for (dim = 0; dim < GridRank; dim++)
-      for (i = 0; i < NumberOfParticles; i++)
-	if (ParticleMass[i] >= 0) {
-	  if (ParticlePosition[dim][i] < GridLeftEdge[dim]) {
-	    NumberToMove[dim*2+0]++;
-	    ParticleMass[i] = -ParticleMass[i];
-	  }
-	  if (ParticlePosition[dim][i] > GridRightEdge[dim]) {
-	    NumberToMove[dim*2+1]++;
-	    ParticleMass[i] = -ParticleMass[i];
-	  }
-	}
- 
-    /* Allocate space. */
- 
-    int TotalToMove = 0, NumberOfParticleFields = 9+NumberOfParticleAttributes;
-    for (i = 0; i < 6; i++) {
-      TotalToMove += NumberToMove[i];
-      if (NumberToMove[i] > 0) {
-//	fprintf(stdout, "P%"ISYM", grid %x: side %"ISYM", NumberToMove = %"ISYM"\n",
-//	       ProcessorNumber, this, i, NumberToMove[i]);
-	fflush(stdout);
-	ParticleData[i] = new float_int[NumberToMove[i]*
-				        NumberOfParticleFields];
-      } else
-	ParticleData[i] = NULL;
-    }
- 
-    /* Set ToGrid. */
- 
-    for (dim = 0; dim < GridRank; dim++) {
-      int DimSize = nint((DomainRightEdge[dim] -
-		      DomainLeftEdge[dim])/CellWidth[dim][0]);
- 
-      /* Find Left grid */
- 
-      for (grid = 0; grid < NumberOfGrids; grid++) {
-	int FoundIt = TRUE;
-	for (i = 0; i < GridRank; i++) {
-	  if (dim != i && nint(GridLeftEdge[i]/CellWidth[i][0]) !=
-	      nint(Grids[grid]->GridLeftEdge[i]/CellWidth[i][0]))
-	    FoundIt = FALSE;
-	  if (dim == i && (nint(
-               Grids[grid]->GridRightEdge[i]/CellWidth[i][0]) % DimSize)
-	      != nint(GridLeftEdge[i]/CellWidth[i][0]))
-	    FoundIt = FALSE;
-	}
-	if (FoundIt) {
-	  ToGrid[dim*2+0] = grid;
-	  break;
-	}
-      }
- 
-      /* Find right grid */
- 
-      for (grid = 0; grid < NumberOfGrids; grid++) {
-	int FoundIt = TRUE;
-	for (i = 0; i < GridRank; i++) {
-	  if (dim != i && nint(GridLeftEdge[i]/CellWidth[i][0]) !=
-	      nint(Grids[grid]->GridLeftEdge[i]/CellWidth[i][0]))
-	    FoundIt = FALSE;
-	  if (dim == i && (nint(
-               GridRightEdge[i]/CellWidth[i][0]) % DimSize)
-	      != nint(Grids[grid]->GridLeftEdge[i]/CellWidth[i][0]))
-	    FoundIt = FALSE;
-	}
-	if (FoundIt) {
-	  ToGrid[dim*2+1] = grid;
-	  break;
-	}
-      }
- 
-    } // end loop over dims
- 
-    if (TotalToMove == 0)
-      return SUCCESS;
- 
-    /* Move particles (mark those moved by setting mass = FLOAT_UNDEFINED). */
- 
-    for (dim = 0; dim < GridRank; dim++) {
-      int n1 = 0, n2 = 0;
-      for (i = 0; i < NumberOfParticles; i++)
-	if (ParticleMass[i] != FLOAT_UNDEFINED) {
- 
-	  /* shift left. */
- 
-	  if (ParticlePosition[dim][i] < GridLeftEdge[dim]) {
-	    for (j = 0; j < GridRank; j++) {
-	      ParticleData[dim*2][n1++].FVAL = ParticlePosition[j][i];
-	      ParticleData[dim*2][n1++].fval = ParticleVelocity[j][i];
-	    }
-	    ParticleData[dim*2][n1++].fval = -ParticleMass[i];
-	    ParticleData[dim*2][n1++].IVAL = ParticleNumber[i];
-	    ParticleData[dim*2][n1++].ival = ParticleType[i];
-	    for (j = 0; j < NumberOfParticleAttributes; j++)
-	      ParticleData[dim*2][n1++].fval = ParticleAttribute[j][i];
-	    ParticleMass[i] = FLOAT_UNDEFINED;
-	  }
- 
-	  /* shift right. */
- 
-	  if (ParticlePosition[dim][i] > GridRightEdge[dim]) {
-	    for (j = 0; j < MAX_DIMENSION; j++) {
-	      ParticleData[dim*2+1][n2++].FVAL = ParticlePosition[j][i];
-	      ParticleData[dim*2+1][n2++].fval = ParticleVelocity[j][i];
-	    }
-	    ParticleData[dim*2+1][n2++].fval = -ParticleMass[i];
-	    ParticleData[dim*2+1][n2++].IVAL = ParticleNumber[i];
-	    ParticleData[dim*2+1][n2++].ival = ParticleType[i];
-	    for (j = 0; j < NumberOfParticleAttributes; j++)
-	      ParticleData[dim*2+1][n2++].fval = ParticleAttribute[j][i];
-	    ParticleMass[i] = FLOAT_UNDEFINED;
-	  }
- 
-	} // end: if (ParticleMass[i] != FLOAT_UNDEFINED)
-    } // end: loop over dims
- 
-  } // end: if (COPY_OUT)
- 
-  /* ----------------------------------------------------------------- */
-  /* Copy particle back into grid. */
- 
-  else {
-
-    double t00, t01, ttt;
-    double T1, T2, T3, T4;
-
-    t00 = 0.0;
-    t01 = 0.0;
-    ttt = 0.0;
-
-    T1 = 0.0;
-    T2 = 0.0;
-    T3 = 0.0;
-    T4 = 0.0;
-
-#ifdef USE_MPI
-    t00 = MPI_Wtime();
-#endif
- 
-    /* Count up total number. */
- 
-    int TotalNumberOfParticles = NumberOfParticles;
- 
-    for (i = 0; i < NumberOfParticles; i++)
-      if (ParticleMass[i] == FLOAT_UNDEFINED)
-	TotalNumberOfParticles--;
- 
-    for (i = 0; i < 6; i++)
-      TotalNumberOfParticles += NumberToMove[i];
- 
-    if (TotalNumberOfParticles == 0 && NumberOfParticles == 0)
-      return SUCCESS;
- 
-    /* Allocate space for the particles. */
-
-#ifdef USE_MPI
-    T1 = MPI_Wtime();
-#endif
- 
-    FLOAT *Position[MAX_DIMENSION];
-    float *Velocity[MAX_DIMENSION], *Mass,
-          *Attribute[MAX_NUMBER_OF_PARTICLE_ATTRIBUTES];
-    PINT *Number;
-    int *Type;
- 
-    Mass = new float[TotalNumberOfParticles];
-    Number = new PINT[TotalNumberOfParticles];
-    Type = new int[TotalNumberOfParticles];
-    for (dim = 0; dim < GridRank; dim++) {
-      Position[dim] = new FLOAT[TotalNumberOfParticles];
-      Velocity[dim] = new float[TotalNumberOfParticles];
-    }
-    for (i = 0; i < NumberOfParticleAttributes; i++)
-      Attribute[i] = new float[TotalNumberOfParticles];
-
-    if (Velocity[GridRank-1] == NULL && TotalNumberOfParticles != 0) {
-      ENZO_FAIL("malloc error (out of memory?)\n");
-    }
-
-#ifdef USE_MPI
-    T2 = MPI_Wtime();
-#endif
- 
-    /* Copy this grid's particles to the new space (don't copy any erased
-       particles, those with Mass == FLOAT_UNDEFINED). */
-
-/* 
-    int n = 0;
-    for (i = 0; i < NumberOfParticles; i++)
-      if (ParticleMass[i] != FLOAT_UNDEFINED) {
-	Mass[n]   = ParticleMass[i];
-	Number[n] = ParticleNumber[i];
-	Type[n] = ParticleType[i];
- 
-	for (dim = 0; dim < GridRank; dim++) {
-	  Position[dim][n] = ParticlePosition[dim][i];
-	  Velocity[dim][n] = ParticleVelocity[dim][i];
-	}
-	for (j = 0; j < NumberOfParticleAttributes; j++)
-	  Attribute[j][n] = ParticleAttribute[j][i];
-	n++;
-      }
-*/
-
-    int n;
-
-    n = 0;
-    for (i = 0; i < NumberOfParticles; i++) {
-      if (ParticleMass[i] != FLOAT_UNDEFINED) {
-        Mass[n]   = ParticleMass[i];
-        Number[n] = ParticleNumber[i];
-        Type[n] = ParticleType[i];
-        n++;
-      }
-    }
-
-    for (dim = 0; dim < GridRank; dim++) {
-      n = 0;
-      for (i = 0; i < NumberOfParticles; i++) {
-        if (ParticleMass[i] != FLOAT_UNDEFINED) {
-          Position[dim][n] = ParticlePosition[dim][i];
-          Velocity[dim][n] = ParticleVelocity[dim][i];
-          n++;
-        }
-      }
-    }
-
-    for (j = 0; j < NumberOfParticleAttributes; j++) {
-      n = 0;
-      for (i = 0; i < NumberOfParticles; i++) {
-        if (ParticleMass[i] != FLOAT_UNDEFINED) {
-          Attribute[j][n] = ParticleAttribute[j][i];
-          n++;
-        }
-      }
-    }
-
-#ifdef USE_MPI
-    T3 = MPI_Wtime();
-#endif
- 
-    /* Copy new particles (and wrap around edges). */
- 
-    for (j = 0; j < 6; j++)
-      if (NumberToMove[j] > 0) {
-	int n2 = 0;
-	for (i = 0; i < NumberToMove[j]; i++) {
- 
-	  for (dim = 0; dim < GridRank; dim++) {
-	    Position[dim][n] = ParticleData[j][n2++].FVAL;
-	    if (Position[dim][n] > DomainRightEdge[dim])
-	      Position[dim][n] -= DomainRightEdge[dim] - DomainLeftEdge[dim];
-	    if (Position[dim][n] < DomainLeftEdge[dim])
-	      Position[dim][n] += DomainRightEdge[dim] - DomainLeftEdge[dim];
-	    Velocity[dim][n] = ParticleData[j][n2++].fval;
-	  }
- 
-	  Mass[n]   = ParticleData[j][n2++].fval;
-	  Number[n] = ParticleData[j][n2++].IVAL;
-	  Type[n] = ParticleData[j][n2++].ival;
-	  for (k = 0; k < NumberOfParticleAttributes; k++)
-	    Attribute[k][n] = ParticleData[j][n2++].fval;
- 
-	  n++;
-	}
-      }
-
-#ifdef USE_MPI
-    T4 = MPI_Wtime();
-#endif
- 
-    /* Set new number of particles in this grid. */
- 
-    NumberOfParticles = TotalNumberOfParticles;
- 
-    /* Delete old pointers and copy new ones into place. */
- 
-    this->DeleteParticles();
-    this->SetParticlePointers(Mass, Number, Type, Position, Velocity,
-			      Attribute);
- 
-#ifdef USE_MPI
-    t01 = MPI_Wtime();
-#endif
-    ttt = t01-t00;
-    // fprintf(stderr, "COPY IN %"ISYM" : %16.6e : %16.6e %16.6e %16.6e %16.6e %16.6e\n", MyProcessorNumber, ttt,
-    //         T1-t00, T2-T1, T3-T2, T4-T3, t01-T4);
-
-  } // end: if (COPY_IN)
-
-
- 
-  return SUCCESS;
-}
-#endif /* OPTIMIZED_CTP */

src/enzo/Grid_CommunicationTransferParticlesOpt.C

 /  PURPOSE:
 /
 ************************************************************************/
-#ifdef OPTIMIZED_CTP
- 
 #ifdef USE_MPI
 #include <mpi.h>
 #endif
 #include "ExternalBoundary.h"
 #include "Grid.h"
 #include "Hierarchy.h"
- 
+
+int search_lower_bound(int *arr, int value, int low, int high, 
+		       int total);
  
 int grid::CommunicationTransferParticles(grid* Grids[], int NumberOfGrids,
-					 int ThisGridNum, int *&NumberToMove, 
+					 int ThisGridNum, int TopGridDims[],
+					 int *&NumberToMove, 
 					 int StartIndex, int EndIndex, 
 					 particle_data *&List,
-					 int *Layout, int *GridMap, 
+					 int *Layout, int *GStartIndex[], 
+					 int *GridMap, 
 					 int CopyDirection)
 {
  
   /* Declarations. */
  
-  int i, j, k, dim, grid, proc, grid_num;
+  int i, j, k, dim, grid, proc, grid_num, bin, CenterIndex;
   int GridPosition[MAX_DIMENSION];
-  int *ToGrid;
+  int *ToGrid, *pbin;
+  FLOAT r[MAX_DIMENSION];
  
   for (dim = 0; dim < MAX_DIMENSION; dim++)
     GridPosition[dim] = 0;
     ToGrid = new int[NumberOfParticles];
 
     float DomainWidth[MAX_DIMENSION], DomainWidthInv[MAX_DIMENSION];
-    float LayoutFloat[MAX_DIMENSION];
     for (dim = 0; dim < MAX_DIMENSION; dim++) {
       DomainWidth[dim] = DomainRightEdge[dim] - DomainLeftEdge[dim];
       DomainWidthInv[dim] = 1.0/DomainWidth[dim];
-      LayoutFloat[dim] = (float) Layout[dim];
     }
 
     // Periodic boundaries
     for (i = 0; i < NumberOfParticles; i++) {
 
       for (dim = 0; dim < GridRank; dim++) {
-	GridPosition[dim] = 
-	  (int) (LayoutFloat[dim] * 
+
+	if (Layout[dim] == 1) {
+	  GridPosition[dim] = 0;
+	} else {
+
+	  CenterIndex = 
+	  (int) (TopGridDims[dim] * 
 		 (ParticlePosition[dim][i] - DomainLeftEdge[dim]) *
 		 DomainWidthInv[dim]);
-	GridPosition[dim] = min(GridPosition[dim], Layout[dim]-1);
-      }
+
+	  GridPosition[dim] = 
+	    search_lower_bound(GStartIndex[dim], CenterIndex, 0, Layout[dim],
+			       Layout[dim]);
+	  GridPosition[dim] = min(GridPosition[dim], Layout[dim]-1);
+	  
+	} // ENDELSE
+
+      } // ENDFOR dim
 
       grid_num = GridPosition[0] + 
 	Layout[0] * (GridPosition[1] + Layout[1]*GridPosition[2]);
+
       grid = GridMap[grid_num];
       if (grid != ThisGridNum) {
 	proc = Grids[grid]->ReturnProcessorNumber();
 	NumberToMove[proc]++;
+//	printf("grid %d->%d: Pos = %f %f %f\n", ThisGridNum, grid, 
+//	       ParticlePosition[0][i], ParticlePosition[1][i], 
+//	       ParticlePosition[2][i]);
       }
+
       ToGrid[i] = grid;
+
     } // ENDFOR particles
 
     /* Allocate space. */
 	  List[n1].proc = MyProcessorNumber;
 	  ParticleMass[i] = FLOAT_UNDEFINED;
 	  n1++;
-	} // ENDIF different processor
+	} // ENDIF different grid
       } // ENDFOR particles
 
     } // ENDIF TotalToMove > PreviousTotalToMove
 
   } // end: if (COPY_IN)
 
-
- 
   return SUCCESS;
 }
-#endif /* OPTIMIZED_CTP */

src/enzo/Grid_CommunicationTransferStars.C

-/***********************************************************************
-/
-/  GRID CLASS (COPY STARS INTO OR OUT OF GRID)
-/
-/  written by: Greg Bryan
-/  date:       January, 1999
-/  modified1:  Robert Harkness
-/  date:       April, 2006
-/  modified2:  John Wise
-/  date:       July, 2009 (modified grid:CTP)
-/
-/  PURPOSE:
-/
-************************************************************************/
-#ifndef OPTIMIZED_CTP
-//
- 
-#ifdef USE_MPI
-#include <mpi.h>
-#endif
-#include <stdio.h>
-
-#include "ErrorExceptions.h"
-#include "macros_and_parameters.h"
-#include "typedefs.h"
-#include "global_data.h"
-#include "Fluxes.h"
-#include "GridList.h"
-#include "ExternalBoundary.h"
-#include "Grid.h"
-#include "Hierarchy.h"
- 
-Star *PopStar(Star * &Node);
-Star* StarBufferToList(StarBuffer *buffer, int n);
-void InsertStarAfter(Star * &Node, Star * &NewNode);
-void DeleteStarList(Star * &Node);
- 
-int grid::CommunicationTransferStars(grid* Grids[], int NumberOfGrids,
-		 int ToGrid[6], int NumberToMove[6],
-		 StarBuffer *StarData[6], int CopyDirection)
-{
- 
-  /* Declarations. */
- 
-  int i, j, k, dim, grid;
-  Star *cstar, *MoveStar, *StarList[6];
- 
-  /* ----------------------------------------------------------------- */
-  /* Copy star out of grid. */
- 
-  if (CopyDirection == COPY_OUT) {
-
-    /* If there are no stars to move, we're done. */
- 
-    if (NumberOfStars == 0)
-      return SUCCESS;
- 
-    /* Count stars to move (mark already counted by setting mass < 0). */
- 
-    for (dim = 0; dim < GridRank; dim++)
-      for (cstar = Stars; cstar; cstar = cstar->NextStar)
-	if (cstar->Mass >= 0) {
-	  if (cstar->pos[dim] < GridLeftEdge[dim]) {
-	    NumberToMove[dim*2+0]++;
-	    cstar->Mass *= -cstar->Mass;
-	  }
-	  if (cstar->pos[dim] > GridRightEdge[dim]) {
-	    NumberToMove[dim*2+1]++;
-	    cstar->Mass = -cstar->Mass;
-	  }
-	}
- 
-    /* Allocate space. */
- 
-    int TotalToMove = 0;
-    for (i = 0; i < 6; i++) {
-      TotalToMove += NumberToMove[i];
-      StarList[i] = NULL;
-      if (NumberToMove[i] > 0)
-	StarData[i] = new StarBuffer[NumberToMove[i]];
-      else
-	StarData[i] = NULL;
-    }
- 
-    /* Set ToGrid. */
- 
-    for (dim = 0; dim < GridRank; dim++) {
-      int DimSize = nint((DomainRightEdge[dim] -
-		      DomainLeftEdge[dim])/CellWidth[dim][0]);
- 
-      /* Find Left grid */
- 
-      for (grid = 0; grid < NumberOfGrids; grid++) {
-	int FoundIt = TRUE;
-	for (i = 0; i < GridRank; i++) {
-	  if (dim != i && nint(GridLeftEdge[i]/CellWidth[i][0]) !=
-	      nint(Grids[grid]->GridLeftEdge[i]/CellWidth[i][0]))
-	    FoundIt = FALSE;
-	  if (dim == i && (nint(
-               Grids[grid]->GridRightEdge[i]/CellWidth[i][0]) % DimSize)
-	      != nint(GridLeftEdge[i]/CellWidth[i][0]))
-	    FoundIt = FALSE;
-	}
-	if (FoundIt) {
-	  ToGrid[dim*2+0] = grid;
-	  break;
-	}
-      }
- 
-      /* Find right grid */
- 
-      for (grid = 0; grid < NumberOfGrids; grid++) {
-	int FoundIt = TRUE;
-	for (i = 0; i < GridRank; i++) {
-	  if (dim != i && nint(GridLeftEdge[i]/CellWidth[i][0]) !=
-	      nint(Grids[grid]->GridLeftEdge[i]/CellWidth[i][0]))
-	    FoundIt = FALSE;
-	  if (dim == i && (nint(
-               GridRightEdge[i]/CellWidth[i][0]) % DimSize)
-	      != nint(Grids[grid]->GridLeftEdge[i]/CellWidth[i][0]))
-	    FoundIt = FALSE;
-	}
-	if (FoundIt) {
-	  ToGrid[dim*2+1] = grid;
-	  break;
-	}
-      }
- 
-    } // end loop over dims
- 
-    if (TotalToMove == 0)
-      return SUCCESS;
- 
-    /* Move stars into linked lists for each face */
-
-    for (dim = 0; dim < GridRank; dim++) {
-      int n1 = 0, n2 = 0;
-      NumberOfStars = 0;
-      cstar = Stars;
-      Stars = NULL;
-      while (cstar != NULL) {
-
-	MoveStar = PopStar(cstar);  // also advances to NextStar
- 
-	/* shift left. */
- 
-	if (MoveStar->pos[dim] < GridLeftEdge[dim]) {
-	  MoveStar->Mass = -MoveStar->Mass;
-	  InsertStarAfter(StarList[dim*2], MoveStar);
-	}
- 
-	/* shift right. */
- 
-	else if (MoveStar->pos[dim] > GridRightEdge[dim]) {
-	  MoveStar->Mass = -MoveStar->Mass;
-	  InsertStarAfter(StarList[dim*2+1], MoveStar);
-	}
-
-	/* no shift */
-
-	else {
-	  InsertStarAfter(this->Stars, MoveStar);
-	  NumberOfStars++;
-	}
- 
-      } // ENDFOR stars
-    } // end: loop over dims
- 
-    /* Convert linked lists into arrays and delete them */
-
-    for (i = 0; i < 6; i++)
-      if (NumberToMove[i] > 0) {
-	StarData[i] = StarList[i]->StarListToBuffer(NumberToMove[i]);
-	DeleteStarList(StarList[i]);
-      }
-
-  } // end: if (COPY_OUT)
- 
-  /* ----------------------------------------------------------------- */
-  /* Copy stars back into grid. */
- 
-  else {
-
-    /* Count up total number. */
- 
-    int TotalNumberOfStars = NumberOfStars;
- 
-    for (i = 0; i < 6; i++)
-      TotalNumberOfStars += NumberToMove[i];
- 
-    if (TotalNumberOfStars == 0 && NumberOfStars == 0)
-      return SUCCESS;
-
-    /* Copy stars from buffer into linked list */
-
-    for (i = 0; i < 6; i++)
-      if (NumberToMove[i] > 0) {
-	MoveStar = StarBufferToList(StarData[i], NumberToMove[i]);
-	MoveStar->CurrentGrid = this;
-	MoveStar->GridID = this->ID;
-	InsertStarAfter(this->Stars, MoveStar);
-      } // ENDIF NumberToMove > 0
- 
-    /* Set new number of stars in this grid. */
- 
-    NumberOfStars = TotalNumberOfStars;
- 
-  } // end: if (COPY_IN)
-
- 
-  return SUCCESS;
-}
-#endif /* OPTIMIZED_CTP */

src/enzo/Grid_CommunicationTransferStarsOpt.C

 /  PURPOSE:
 /
 ************************************************************************/
-#ifdef OPTIMIZED_CTP
- 
 #ifdef USE_MPI
 #include <mpi.h>
 #endif
 Star *PopStar(Star * &Node);
 Star* StarBufferToList(StarBuffer buffer);
 void InsertStarAfter(Star * &Node, Star * &NewNode);
+int search_lower_bound(int *arr, int value, int low, int high, 
+		       int total);
  
 int grid::CommunicationTransferStars(grid* Grids[], int NumberOfGrids,
-				     int ThisGridNum, int *&NumberToMove, 
+				     int ThisGridNum, int TopGridDims[],
+				     int *&NumberToMove, 
 				     int StartIndex, int EndIndex, 
-				     star_data *&List,
-				     int *Layout, int *GridMap, 
+				     star_data *&List, int *Layout, 
+				     int *GStartIndex[], int *GridMap, 
 				     int CopyDirection)
 {
  
   /* Declarations. */
  
-  int i, j, k, dim, grid, proc, grid_num;
+  int i, j, k, dim, grid, proc, grid_num, width, bin, CenterIndex;
   int GridPosition[MAX_DIMENSION];
-  int *ToGrid;
+  FLOAT r[MAX_DIMENSION];
+  int *ToGrid, *pbin;
   Star *cstar, *MoveStar;
 
   for (dim = 0; dim < MAX_DIMENSION; dim++)
     ToGrid = new int[NumberOfStars];
 
     float DomainWidth[MAX_DIMENSION], DomainWidthInv[MAX_DIMENSION];
-    float LayoutFloat[MAX_DIMENSION];
     for (dim = 0; dim < MAX_DIMENSION; dim++) {
       DomainWidth[dim] = DomainRightEdge[dim] - DomainLeftEdge[dim];
       DomainWidthInv[dim] = 1.0/DomainWidth[dim];
-      LayoutFloat[dim] = (float) Layout[dim];
     }
 
     // Periodic boundaries
     for (cstar = Stars, i = 0; cstar; cstar = cstar->NextStar, i++) {
 
       for (dim = 0; dim < GridRank; dim++) {