Permalink
Browse files

Fixed bug where OUTPUT_POTENTIAL conflicted with SINK_PARTICLES. Vari…

…ous printf statements currently poluting the code from the debugging.
  • Loading branch information...
1 parent 3fb3ffd commit 103db8ab660e0d61764316d2f6b3da811eee3b2e Matthew Young committed May 28, 2012
@@ -302,6 +302,7 @@ extern struct global_data_all_processes
#ifdef SINK_PARTICLES
/* accretion and wind parameters */
double AccretionRadius; /*!< the radius at which sph particles are accreted onto sinks */
+ int AccreteFlag; /* Determines whether a domain decomposition should cause accretion to be processed */
#endif
/* some force counters */
@@ -114,12 +114,12 @@ void domain_Decomposition(void)
#ifdef SINK_PARTICLES
- if(All.NumCurrentTiStep > 0 ){
+ if(All.NumCurrentTiStep > 0 && All.AccreteFlag){
// first see if any processors have particles scheduled for accretion
MPI_Barrier(MPI_COMM_WORLD);
MPI_Allreduce(&AccNum, &AccNumTot, 1, MPI_INT, MPI_SUM, MPI_COMM_WORLD);
- // if so, accrete the,
+ // if so, accrete them
if(AccNumTot > 0){
destroy_doomed_particles();
MPI_Barrier(MPI_COMM_WORLD);
@@ -28,7 +28,7 @@ void gravity_tree(void)
{
long long ntot;
int numnodes, nexportsum = 0;
- int i, j, iter = 0;
+ int i, j, iter = 0,ndonetot=0,nloop=0;
int *numnodeslist, maxnumnodes, nexport, *numlist, *nrecv, *ndonelist;
double tstart, tend, timetree = 0, timecommsumm = 0, timeimbalance = 0, sumimbalance;
double ewaldcount;
@@ -99,20 +99,53 @@ void gravity_tree(void)
i = 0; /* beginn with this index */
ntotleft = ntot; /* particles left for all tasks together */
+ ndonetot = 0;
+ MPI_Allgather(&NumPart, 1, MPI_INT, ndonelist, 1, MPI_INT, MPI_COMM_WORLD);
+ for(j = 0; j < NTask; j++)
+ ndonetot += ndonelist[j];
+ if(ThisTask==0)
+ printf("The total number of particles is %d\n",ndonetot);
+
+
+ if(ThisTask ==0)
+ {
+ printf("Beginning non-trivial stuff.\n");
+ printf("There are %d particles marked for accretion.\n",AccNum);
+ printf("There are %d particles that need processing.\n",ntot);
+ }
+ for(j=0;j<NumPart;j++)
+ {
+ if(P[j].ID==4678)
+ {
+ printf("At time %d, bad particle has (Start,End) = (%d,%d), Acc = %d\n",All.Ti_Current,P[j].Ti_begstep,P[j].Ti_endstep,SphP[j].AccretionTarget);
+ }
+ }
while(ntotleft > 0)
{
+ if(ThisTask ==0)
+ printf("Begin the while loop. With %d particle to do. NumPart is %d. Starting at index %d\n",ntotleft,NumPart,i);
iter++;
for(j = 0; j < NTask; j++)
nsend_local[j] = 0;
+ for(j=0;j<NumPart;j++)
+ {
+ if(P[j].ID==4678)
+ {
+ printf("At time %d, before main loop, particle has (Start,End) = (%d,%d), Acc = %d\n",All.Ti_Current,P[j].Ti_begstep,P[j].Ti_endstep,SphP[j].AccretionTarget);
+ }
+ }
+
/* do local particles and prepare export list */
tstart = second();
for(nexport = 0, ndone = 0; i < NumPart && nexport < All.BunchSizeForce - NTask; i++)
+ {
+ //printf("Time is %d\n",All.Ti_Current);
+ nloop++;
if(P[i].Ti_endstep == All.Ti_Current)
{
ndone++;
-
for(j = 0; j < NTask; j++)
Exportflag[j] = 0;
#ifndef PMGRID
@@ -149,7 +182,26 @@ void gravity_tree(void)
nsend_local[j]++;
}
}
- }
+ }else{
+ printf("In the main loop, it has ID = %d, Type %d, (start,end) = (%d,%d) and current = %d\n",P[i].ID,P[i].Type,P[i].Ti_begstep,P[i].Ti_endstep,All.Ti_Current);
+ }
+
+ }
+ ndonetot = 0;
+ MPI_Allgather(&ndone, 1, MPI_INT, ndonelist, 1, MPI_INT, MPI_COMM_WORLD);
+ for(j = 0; j < NTask; j++)
+ ndonetot += ndonelist[j];
+ if(ThisTask==0)
+ printf("We have processed %d particles.\n",ndonetot);
+ ndonetot = 0;
+ MPI_Allgather(&nloop, 1, MPI_INT, ndonelist, 1, MPI_INT, MPI_COMM_WORLD);
+ for(j = 0; j < NTask; j++)
+ ndonetot += ndonelist[j];
+ if(ThisTask==0)
+ printf("But the loop has seen %d particles.\n",ndonetot);
+
+ if(ThisTask==0)
+ printf("Local particles done. Processed %d particles. Need to export %d of them.\n",ndone,nexport);
tend = second();
timetree += timediff(tstart, tend);
@@ -223,7 +275,9 @@ void gravity_tree(void)
}
tend = second();
timetree += timediff(tstart, tend);
-
+ if(ThisTask==0)
+ printf("Remote particles done.\n");
+
tstart = second();
MPI_Barrier(MPI_COMM_WORLD);
tend = second();
@@ -281,11 +335,15 @@ void gravity_tree(void)
level = ngrp - 1;
}
-
+ if(ThisTask==0)
+ printf("particles done.\n");
+
MPI_Allgather(&ndone, 1, MPI_INT, ndonelist, 1, MPI_INT, MPI_COMM_WORLD);
for(j = 0; j < NTask; j++)
ntotleft -= ndonelist[j];
}
+ if(ThisTask==0)
+ printf("Done all but trivial tasks. ntotleft is %d\n",ntotleft);
free(ndonelist);
free(nsend);
@@ -23,6 +23,9 @@ void init(void)
double a3;
All.Time = All.TimeBegin;
+#ifdef SINK_PARTICLES
+ All.AccreteFlag = 1;
+#endif
switch (All.ICFormat)
{
@@ -165,7 +165,7 @@ void identify_doomed_particles(void)
double *local_sink_velx, *local_sink_vely, *local_sink_velz;
double *local_sink_mass;
int *local_sink_ID;
- int verbose = 0;
+ int verbose = 1;
FLOAT *pos, *vel;
FLOAT Postemp[3], Veltemp[3];
@@ -483,7 +483,6 @@ void destroy_doomed_particles(void)
P[j].ID,ThisTask,AccNum,P[j].Vel[0],P[j].Vel[1],P[j].Vel[2], \
P[j].Pos[0],P[j].Pos[1],P[j].Pos[2],P[j].Mass);
P[j].Ti_endstep = All.Ti_Current;
-
}
}
}
@@ -500,6 +499,8 @@ void destroy_doomed_particles(void)
if(P[i].Ti_endstep == All.Ti_Current){
NumForceUpdate--;
NumSphUpdate--;
+ }else{
+ printf("This is not a current time step particle.\nIt has ID,Type,start,end,current %d,%d,%d,%d,%d\n",P[i].ID,P[i].Type,P[i].Ti_begstep,P[i].Ti_endstep,All.Ti_Current);
}
for(k = i+1; k<=NumPart; k++){ // actually remove the particle here,
// and shift everything down to fill the gap in the array.
@@ -60,6 +60,7 @@ void run(void)
identify_doomed_particles();
MPI_Allreduce(&AccNum, &AccNumTot, 1, MPI_INT, MPI_SUM, MPI_COMM_WORLD);
+ //This will cause domain decomposition to be performed next time this loop is iterated
if(AccNumTot > 0) All.NumForcesSinceLastDomainDecomp = All.TotNumPart * All.TreeDomainUpdateFrequency + 1;
}
#endif
@@ -232,7 +233,10 @@ void find_next_sync_point_and_drift(void)
#ifdef OUTPUTPOTENTIAL
All.NumForcesSinceLastDomainDecomp = 1 + All.TotNumPart * All.TreeDomainUpdateFrequency;
+ //Don't want to accrete here, want to do that when we expect it.
+ All.AccreteFlag=0;
domain_Decomposition();
+ All.AccreteFlag=1;
compute_potential();
#endif
savepositions(All.SnapshotFileCount++); /* write snapshot file */
Binary file not shown.
@@ -1,10 +1,10 @@
%Filenames and file formats
OutputDir /data/my304/Output/Tests
-SnapshotFileBase ThinDisc
+SnapshotFileBase LowResThinDisc
SnapFormat 1
NumFilesPerSnapshot 1
-InitCondFile /data/my304/Working/gadget2-AccretionDiscs/Gadget2AccretionDiscs/Tests/thinDisc.IC
+InitCondFile /data/my304/Working/gadget2-AccretionDiscs/Gadget2AccretionDiscs/Tests/LowResThinDisc.IC
ICFormat 1
RestartFile restart
InfoFile info.txt
@@ -60,8 +60,8 @@ MaxRMSDisplacementFac 0.25
OutputListOn 0
OutputListFilename output_times.txt
TimeOfFirstSnapshot 0.0
-TimeBetSnapshot 1
-TimeBetStatistics 1
+TimeBetSnapshot .025
+TimeBetStatistics .025
NumFilesWrittenInParallel 1
%Unit system

0 comments on commit 103db8a

Please sign in to comment.