Hi,
if I understand correctly you are running the same code twice to get the
two accelerations? I would recommend to dump both of the accelerations
to the same snapshot, requiring you to run the code only once. The
cleanest way to do this would be to add two new IO blocks, e.g.
IO_GRAVACCEL and IO_HYDROACCEL. You then just have to scan the io.c file
for each IO_ACCEL location and add IO_GRAVACCEL and IO_HYDROACCEL in an
analogous way. Also remember to update allvars.h: add the two blocks +
update the number of IO blocks just above the enum.
As for your current implementation: I must say that I don't see anything
wrong with it. You could try using a debugger or print commands to check
what is actually happening inside the code, e.g. add:
if (P[pindex].ID == 0)
{
printf("(%g,%g,%g) (%g,%g,%g) (%g,%g,%g) (%g,%g,%g)\n",
fac1*P[pindex].GravAccel[0],fac1*P[pindex].GravAccel[1],fac1*P[pindex].GravAccel[2],
fac1 * P[pindex].GravPM[0],fac1 * P[pindex].GravPM[1],fac1 *
P[pindex].GravPM[2],
fac2*SphP[pindex].HydroAccel[0],fac2*SphP[pindex].HydroAccel[1],fac2*SphP[pindex].HydroAccel[2],
fp[0], fp[1], fp[2]);
}
If the results from that are what you would expect (e.g. fp is different
in the 2 codes) then maybe the problem is located in the code reading
the snapshots.
cheerz,
Sander
Guillermo Arreaga Garcia wrote:
> Hi all,
>
> I need to obtain separately the array of hydrodynamical and gravitational
> accelerations of a gadget2īs simulation having only gas particles. The
> procedure what I have followed is outlined below.
>
> Firstly, I have already activated the option -DOUTPUTACCELERATION in the
> Makefile in order to write the 3-array P[i].Acce[k] for every
> particle i to the output snapshot file.
>
> Secondly, I have modified the file io.c only in the part where
> it will write the particle data to the buffer. My modification consisted in
>
> multiplying by 0.0 the part of the acceleration which I would
> like to eliminate from the array of total acceleration.
>
> That is:
>
> a) to get hydrodinamical acceleration
>
> case IO_ACCEL:
> #ifdef OUTPUTACCELERATION
> for(n = 0; n < pc; pindex++)
> if(P[pindex].Type == type)
> {
> for(k = 0; k < 3; k++)
> fp[k] = 0.0*fac1*P[pindex].GravAccel[k];
> #ifdef PMGRID
> for(k = 0; k < 3; k++)
> fp[k] += 0.0*fac1 * P[pindex].GravPM[k];
> #endif
> if(P[pindex].Type == 0)
> for(k = 0; k < 3; k++)
> fp[k] += fac2*SphP[pindex].HydroAccel[k];
> fp += 3;
> n++;
> }
> #endif
> break;
>
> and
>
> b) to get gravitational acceleration
>
> case IO_ACCEL: /* acceleration */
> #ifdef OUTPUTACCELERATION
> for(n = 0; n < pc; pindex++)
> if(P[pindex].Type == type)
> {
> for(k = 0; k < 3; k++)
> fp[k] = fac1 * P[pindex].GravAccel[k];
> #ifdef PMGRID
> for(k = 0; k < 3; k++)
> fp[k] += fac1 * P[pindex].GravPM[k];
> #endif
> if(P[pindex].Type == 0)
> for(k = 0; k < 3; k++)
> fp[k] += 0.0*fac2*SphP[pindex].HydroAccel[k];
> fp += 3;
> n++;
> }
> #endif
> break;
>
> Thirdly, for the two cases a and b, I update the file io.c with the
> appropriate changes and then I compile and run the gadget2 separately. I have
> written a code for reading the output snapshot file from which I get the Acce
> array.
>
> However, I am surprising because what I get in the two cases are
> the same data for the Acce array...that is, it seems that there is
> no difference between the hydrodynamical and gravitational
> accelerations, which is impossible.
>
> Then I wonder what I am missing or doing incorrectly in the
> previous procedure.
>
> I would greatly appreciate any comment or assistance.
>
>
> Sincerely,
>
> Guillermo Arreaga,
>
>
>
>
> -----------------------------------------------------------
>
> 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
>
--
--------------------------------------------------
| Sander Valcke |
| Sterrenkundig Observatorium, Universiteit Gent |
| Krijgslaan 281 S9, B-9000 Gent, Belgium |
| phone +32-9-264-4796 | fax +32-9-264-4989 |
| https://sourceforge.net/projects/hyplot/ |
--------------------------------------------------
Sometimes the best medicine is to stop taking something.
--------------------------------------------------
Received on 2009-08-27 10:56:55