Hi Eugene,
I think the reason of this problem is that some of your particles have
left the region that is covered by the mesh. In the force computation,
gadget2 recovers from this problem when it occurs by redetermining the
mesh extent and by recomputing the non-periodic kernel. However, the
potential energy computation (which may occur at times in between two
force calculations) did not have this auto-recovery so far.
One possible work-around is use a reasonable value for ENLARGEREGION
(and not 1.0 like you did, but rather 1.2 or so). However, this is not
completely safe, especially when your particle region expands quite a bit.
Another possibility is to add the `auto-readjust' also to the
PM-potential computation. For this, three things need to be changed in
the code:
(1) In potential.c, the lines
#ifdef PMGRID
#ifdef PERIODIC
pmpotential_periodic();
#ifdef PLACEHIGHRESREGION
pmpotential_nonperiodic(1);
#endif
#else
pmpotential_nonperiodic(0);
#ifdef PLACEHIGHRESREGION
pmpotential_nonperiodic(1);
#endif
#endif
#endif
should be replaced by something like this:
#ifdef PMGRID
#ifdef PERIODIC
pmpotential_periodic();
#ifdef PLACEHIGHRESREGION
i = pmpotential_nonperiodic(1);
if(i == 1)
{
pm_init_regionsize();
pm_setup_nonperiodic_kernel();
i = pmpotential_nonperiodic(1);
}
if(i == 1)
endrun(78686);
#endif
#else
pmpotential_nonperiodic(0);
#ifdef PLACEHIGHRESREGION
i = pmpotential_nonperiodic(1);
if(i == 1)
{
pm_init_regionsize();
pm_setup_nonperiodic_kernel();
i = pmpotential_nonperiodic(1);
}
if(i == 1)
endrun(78687);
#endif
#endif
#endif
(2) The return type of the function pmpotential_nonperiodic() needs to
be changed in proto.h from void to int
(3) In pm_nonperiodic.c, the function pmpotential_nonperiodic() needs to
be changed a bit at the beginning. The return type needs to be changed
to int, and two local int-variables, 'flag' and 'flagsum', should be
defined. In the body of the first loop "for(i = 0; i < NumPart; i++)"
one should add the statements
#ifdef PLACEHIGHRESREGION
if(grnr == 0 || (grnr == 1 && ((1 << P[i].Type) &
(PLACEHIGHRESREGION))))
#endif
{
for(j = 0; j < 3; j++)
{
if(P[i].Pos[j] < All.Xmintot[grnr][j] || P[i].Pos[j] >
All.Xmaxtot[grnr][j])
{
if(flag == 0)
{
printf
("Particle Id=%d on task=%d with coordinates (%g|%g|%g) lies outside
PM mesh.\n",
(int)P[i].ID, ThisTask, P[i].Pos[0], P[i].Pos[1], P[i].Pos[2]);
fflush(stdout);
}
flag++;
break;
}
}
}
After the closing bracket of the loop, one should add
MPI_Allreduce(&flag, &flagsum, 1, MPI_INT, MPI_SUM, MPI_COMM_WORLD);
if(flagsum > 0)
{
if(ThisTask == 0)
{
printf("In total %d particles were outside allowed range.\n", flagsum);
fflush(stdout);
}
return(1);
}
And finally, at the very end of the function, a "return 0;" should be added.
Volker
Eugene Kurbatov wrote:
>
> Hello!
>
> I've compiled Gadget-2 (publicly available version) with such options:
> OPT += -DUNEQUALSOFTENINGS
> OPT += -DPEANOHILBERT
> OPT += -DWALLCLOCK
> OPT += -DPMGRID=128
> OPT += -DSYNCRONIZATION
> OPT += -DOUTPUTPOTENTIAL
> OPT += -DCOMPUTE_POTENTIAL_ENERGY
>
> In a snapshot time the run ends with the exit code 131231. If I unset
> last two options then run goes on well. I think, the problem is in the
> `pmpotential_nonperiodic' procedure but I don't know how to defeat it.
>
> It's also possible to set off the PMGRID and expand the ASMTH value (if
> I understood it right) but it is not the way.
>
> Can anybody comment it?
>
> Thanks!
>
>
>
>
> -----------------------------------------------------------
>
> If you wish to unsubscribe from this mailing, send mail to
> minimalist_at_MPA-Garching.MPG.de with a subject of: unsubscribe gadget-list
> A web-archive of this mailing list is available here:
> http://www.mpa-garching.mpg.de/gadget/gadget-list
Received on 2005-11-28 18:02:38