Automated merge with http://hg.spacepope.org/week-of-code

--HG--
branch : week-of-code

src/enzo/CommunicationCollectParticles.C

@@ -53,7 +53,7 @@ int CommunicationShareStars(int *NumberToMove, star_data* &SendList,
 
 int CommunicationCollectParticles(LevelHierarchyEntry *LevelArray[],
 				  int level, bool ParticlesAreLocal, 
-				  int CollectMode)
+				  bool SyncNumberOfParticles, int CollectMode)
 {
 
   /* Create pointer arrays and count grids */
@@ -117,15 +117,19 @@ int CommunicationCollectParticles(LevelHierarchyEntry *LevelArray[],
 
 #ifdef DEBUG_CCP
     for (i = 0; i < NumberOfGrids; i++)
-      printf("CCP[P%"ISYM"B]: grid %"ISYM", %"ISYM" proc, %"ISYM" particles\n",
+      printf("CCP[P%"ISYM"B]: grid %"ISYM", %"ISYM" proc, %"ISYM" particles, "
+	     "%"ISYM" stars\n",
 	     MyProcessorNumber, i,
 	     GridHierarchyPointer[i]->GridData->ReturnProcessorNumber(),
-	     GridHierarchyPointer[i]->GridData->ReturnNumberOfParticles());
+	     GridHierarchyPointer[i]->GridData->ReturnNumberOfParticles(),
+	     GridHierarchyPointer[i]->GridData->ReturnNumberOfStars());
     for (i = 0; i < NumberOfSubgrids; i++)
-      printf("CCP[P%"ISYM"B]: subgrid %"ISYM", %"ISYM" proc, %"ISYM" particles\n",
+      printf("CCP[P%"ISYM"B]: subgrid %"ISYM", %"ISYM" proc, %"ISYM" particles, "
+	     "%"ISYM" stars\n",
 	     MyProcessorNumber, i,
 	     SubgridHierarchyPointer[i]->GridData->ReturnProcessorNumber(),
-	     SubgridHierarchyPointer[i]->GridData->ReturnNumberOfParticles());
+	     SubgridHierarchyPointer[i]->GridData->ReturnNumberOfParticles(),
+	     SubgridHierarchyPointer[i]->GridData->ReturnNumberOfStars());
 #endif /* DEBUG_CCP */
 
     for (i = 0; i < NumberOfProcessors; i++) {
@@ -246,7 +250,8 @@ int CommunicationCollectParticles(LevelHierarchyEntry *LevelArray[],
        processor, set number of particles so everybody agrees. 
     ************************************************************************/
 
-    if ((!KeepLocal && NumberOfProcessors > 1) || ParticlesAreLocal) {
+    if ((!KeepLocal && NumberOfProcessors > 1) || 
+	(ParticlesAreLocal && SyncNumberOfParticles)) {
 
       CommunicationSyncNumberOfParticles(GridHierarchyPointer, NumberOfGrids);
       CommunicationSyncNumberOfParticles(SubgridHierarchyPointer, NumberOfSubgrids);
@@ -255,15 +260,19 @@ int CommunicationCollectParticles(LevelHierarchyEntry *LevelArray[],
 
 #ifdef DEBUG_CCP
     for (i = 0; i < NumberOfGrids; i++)
-      printf("CCP[P%"ISYM"A]: grid %"ISYM", %"ISYM" proc, %"ISYM" particles\n",
+      printf("CCP[P%"ISYM"A]: grid %"ISYM", %"ISYM" proc, %"ISYM" particles, "
+	     "%"ISYM" stars\n",
 	     MyProcessorNumber, i,
 	     GridHierarchyPointer[i]->GridData->ReturnProcessorNumber(),
-	     GridHierarchyPointer[i]->GridData->ReturnNumberOfParticles());
+	     GridHierarchyPointer[i]->GridData->ReturnNumberOfParticles(),
+	     GridHierarchyPointer[i]->GridData->ReturnNumberOfStars());
     for (i = 0; i < NumberOfSubgrids; i++)
-      printf("CCP[P%"ISYM"A]: subgrid %"ISYM", %"ISYM" proc, %"ISYM" particles\n",
+      printf("CCP[P%"ISYM"A]: subgrid %"ISYM", %"ISYM" proc, %"ISYM" particles, "
+	     "%"ISYM" stars\n",
 	     MyProcessorNumber, i,
 	     SubgridHierarchyPointer[i]->GridData->ReturnProcessorNumber(),
-	     SubgridHierarchyPointer[i]->GridData->ReturnNumberOfParticles());
+	     SubgridHierarchyPointer[i]->GridData->ReturnNumberOfParticles(),
+	     SubgridHierarchyPointer[i]->GridData->ReturnNumberOfStars());
 #endif /* DEBUG_CCP */
 
     /* Cleanup. */
@@ -444,8 +453,17 @@ int CommunicationCollectParticles(LevelHierarchyEntry *LevelArray[],
        If the particles and stars are only on the grid's host
        processor, set number of particles so everybody agrees. 
     ************************************************************************/
-
-    CommunicationSyncNumberOfParticles(GridHierarchyPointer, NumberOfGrids);
+
+    if (SyncNumberOfParticles)
+      CommunicationSyncNumberOfParticles(GridHierarchyPointer, NumberOfGrids);
+    else {
+      for (i = StartGrid; i < EndGrid; i++)
+	if (MyProcessorNumber != 
+	    GridHierarchyPointer[i]->GridData->ReturnProcessorNumber()) {
+	  GridHierarchyPointer[i]->GridData->SetNumberOfParticles(0);
+	  GridHierarchyPointer[i]->GridData->SetNumberOfStars(0);
+	}
+    }
 
 #ifdef DEBUG_CCP
     for (i = StartGrid; i < EndGrid; i++)

src/enzo/CommunicationTransferParticles.C

@@ -344,25 +344,27 @@ int CommunicationTransferParticles(grid *GridPointer[], int NumberOfGrids)
     } // end: if grid is on my processor
 
   /* Set number of particles 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]->SetNumberOfParticles(
-//		         GridPointer[j]->ReturnNumberOfParticles()+Changes[j]);
-//      //      printf("Pb(%"ISYM") CTP grid[%"ISYM"] = %"ISYM"\n", MyProcessorNumber, j, GridPointer[j]->ReturnNumberOfParticles());
-//    }
-//    delete [] Changes;
-//  }
+
+#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. */
 

src/enzo/CommunicationTransferParticlesOpt.C

@@ -31,17 +31,18 @@
 #include "TopGridData.h"
 #include "Hierarchy.h"
 #include "LevelHierarchy.h"
+#include "CommunicationUtilities.h"
 void my_exit(int status);
 
 // function prototypes
 
 Eint32 compare_grid(const void *a, const void *b);
 int Enzo_Dims_create(int nnodes, int ndims, int *dims);
-int CommunicationAllSumIntegerValues(int *Values, int Number);
 int CommunicationShareParticles(int *NumberToMove, particle_data* &SendList,
 				int &NumberOfReceives,
 				particle_data* &SharedList);
 
+#define NO_DEBUG_CTP
 #define KEEP_PARTICLES_LOCAL
 
 int CommunicationTransferParticles(grid *GridPointer[], int NumberOfGrids)
@@ -190,23 +191,26 @@ int CommunicationTransferParticles(grid *GridPointer[], int NumberOfGrids)
 
   /* Set number of particles so everybody agrees. */
 
-  // (JHW, May 2009) No longer needed because we synchronize in
-  // CommunicationCollectParticles now.
-
-//  if (NumberOfProcessors > 1) {
-//    int *AllNumberOfParticles = new int[NumberOfGrids];
-//    for (j = 0; j < NumberOfGrids; j++)
-//      if (GridPointer[j]->ReturnProcessorNumber() == MyProcessorNumber)
-//	AllNumberOfParticles[j] = GridPointer[j]->ReturnNumberOfParticles();
-//      else
-//	AllNumberOfParticles[j] = 0;
-//
-//    CommunicationAllSumIntegerValues(AllNumberOfParticles, NumberOfGrids);
-//    for (j = 0; j < NumberOfGrids; j++)
-//      GridPointer[j]->SetNumberOfParticles(AllNumberOfParticles[j]);
-//
-//    delete [] AllNumberOfParticles;
-//  }
+#ifdef UNUSED
+  if (NumberOfProcessors > 1) {
+    int *AllNumberOfParticles = new int[NumberOfGrids];
+    for (j = 0; j < NumberOfGrids; j++)
+      if (GridPointer[j]->ReturnProcessorNumber() == MyProcessorNumber)
+	AllNumberOfParticles[j] = GridPointer[j]->ReturnNumberOfParticles();
+      else
+	AllNumberOfParticles[j] = 0;
+
+    CommunicationAllSumValues(AllNumberOfParticles, NumberOfGrids);
+    printf("P%d: npart = %d %d (%d %d)\n", MyProcessorNumber,
+	   GridPointer[0]->ReturnNumberOfParticles(),
+	   GridPointer[1]->ReturnNumberOfParticles(),
+	   AllNumberOfParticles[0], AllNumberOfParticles[1]);
+    for (j = 0; j < NumberOfGrids; j++)
+      GridPointer[j]->SetNumberOfParticles(AllNumberOfParticles[j]);
+
+    delete [] AllNumberOfParticles;
+  }
+#endif /* UNUSED */
 
   /* Cleanup. */
 
@@ -216,11 +220,7 @@ int CommunicationTransferParticles(grid *GridPointer[], int NumberOfGrids)
   delete [] NumberToMove;
   delete [] GridMap;
 
-#ifdef USE_MPI
-  int temp_int = TotalNumberToMove;
-  MPI_Reduce(&temp_int, &TotalNumberToMove, 1, IntDataType, MPI_SUM, 
-	     ROOT_PROCESSOR, MPI_COMM_WORLD);
-#endif
+  CommunicationSumValues(&TotalNumberToMove, 1);
   if (debug)
     printf("CommunicationTransferParticles: moved = %"ISYM"\n",
   	   TotalNumberToMove);

src/enzo/CommunicationTransferStars.C

@@ -335,6 +335,29 @@ int CommunicationTransferStars(grid *GridPointer[], int NumberOfGrids)
 	 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 */
 
   /* CleanUp. */
 

src/enzo/CommunicationTransferStarsOpt.C

@@ -32,13 +32,13 @@
 #include "TopGridData.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 CommunicationAllSumIntegerValues(int *Values, int Number);
 
 // Remove define.  This method will always be used.
 //#define KEEP_PARTICLES_LOCAL
@@ -137,6 +137,26 @@ int CommunicationTransferStars(grid *GridPointer[], int NumberOfGrids)
     } // ENDFOR grids
   } // ENDIF NumberOfRecieves > 0
 
+
+  /* Set number of stars so everybody agrees. */
+
+#ifdef UNUSED
+  if (NumberOfProcessors > 1) {
+    int *AllNumberOfStars = new int[NumberOfGrids];
+    for (j = 0; j < NumberOfGrids; j++)
+      if (GridPointer[j]->ReturnProcessorNumber() == MyProcessorNumber)
+	AllNumberOfStars[j] = GridPointer[j]->ReturnNumberOfStars();
+      else
+	AllNumberOfStars[j] = 0;
+
+    CommunicationAllSumValues(AllNumberOfStars, NumberOfGrids);
+    for (j = 0; j < NumberOfGrids; j++)
+      GridPointer[j]->SetNumberOfStars(AllNumberOfStars[j]);
+
+    delete [] AllNumberOfStars;
+  }
+#endif
+
   /* Cleanup. */
 
   if (SendList != SharedList)
@@ -145,11 +165,7 @@ int CommunicationTransferStars(grid *GridPointer[], int NumberOfGrids)
   delete [] NumberToMove;
   delete [] GridMap;
 
-#ifdef USE_MPI
-  int temp_int = TotalNumberToMove;
-  MPI_Reduce(&temp_int, &TotalNumberToMove, 1, IntDataType, MPI_SUM, 
-	     ROOT_PROCESSOR, MPI_COMM_WORLD);
-#endif
+  CommunicationSumValues(&TotalNumberToMove, 1);
   if (debug)
     printf("CommunicationTransferStars: moved = %"ISYM"\n",
   	   TotalNumberToMove);

src/enzo/ExternalBoundary.h

@@ -126,6 +126,8 @@ class ExternalBoundary
   int DetachForcingFromBaryonFields();
 
   int AddField(int FieldType);
+  int DeleteObsoleteFields(int *ObsoleteFields, 
+			   int NumberOfObsoleteFields);
 
 };
 

src/enzo/ExternalBoundary_AddField.C

@@ -1,9 +1,9 @@
 /***********************************************************************
 /
-/  EXTERNAL BOUNDARY CLASS (READ THE EXTERNAL BOUNDARY VALUES)
+/  EXTERNAL BOUNDARY CLASS (ADD EXTERNAL BOUNDARY VALUES FOR A NEW FIELD)
 /
-/  written by: Greg Bryan
-/  date:       November, 1994
+/  written by: John Wise
+/  date:       March, 2009
 /  modified1:
 /
 /  PURPOSE:

src/enzo/Grid.h

@@ -1732,6 +1732,8 @@ int CollapseTestInitializeGrid(int NumberOfSpheres,
   int ReadRandomForcingFields(FILE *main_file_pointer);
 
   int AddFields(int TypesToAdd[], int NumberOfFields);
+  int DeleteObsoleteFields(int *ObsoleteFields, 
+			   int NumberOfObsoleteFields);
 
   inline bool isLocal () {return MyProcessorNumber == ProcessorNumber; };
 
@@ -1968,9 +1970,8 @@ int CollapseTestInitializeGrid(int NumberOfSpheres,
 // Radiative transfer methods that don't fit in the TRANSFER define
 //------------------------------------------------------------------------
 
-  int IdentifyRadiativeTransferFields(int &kphHINum, int &gammaHINum,
-				      int &kphHeINum, int &gammaHeINum,
-				      int &kphHeIINum, int &gammaHeIINum,
+  int IdentifyRadiativeTransferFields(int &kphHINum, int &gammaNum,
+				      int &kphHeINum, int &kphHeIINum, 
 				      int &kdissH2INum);
 
 #ifdef TRANSFER

src/enzo/Grid_AddH2Dissociation.C

@@ -60,13 +60,9 @@ int grid::AddH2Dissociation(Star *AllStars)
   /* Get photo-ionization fields */
 
   int kphHINum, kphHeINum, kphHeIINum, kdissH2INum;
-  int gammaHINum, gammaHeINum, gammaHeIINum;
-  if (IdentifyRadiativeTransferFields(kphHINum, gammaHINum, kphHeINum, 
-				      gammaHeINum, kphHeIINum, gammaHeIINum, 
-				      kdissH2INum) == FAIL) {
-    fprintf(stderr, "Error in grid->IdentifyRadiativeTransferFields.\n");
-    ENZO_FAIL("");
-  }
+  int gammaNum;
+  IdentifyRadiativeTransferFields(kphHINum, gammaNum, kphHeINum, kphHeIINum, 
+				  kdissH2INum);
 
   /* For now, initialize H2 photo-dissociation field. */