Re: Extracting the potential on a mesh

From: Volker Springel <volker_at_MPA-Garching.MPG.DE>
Date: Thu, 15 Sep 2005 16:50:08 +0200

Wolfgang Kapferer wrote:
> Dear all,
>
> has someone an easy and fast way to get a file, containing the Potential
> on the mesh defined by PMGRID, calculated by pmpotential_periodic in
> pm_periodic.c. Would MPI_GATHER for all rhogrids in each node be
> sufficient? Is there more than the potential in the arrays rhogrid stored?
>
> Thanks for any useful hints.
>
> Wolfgang
>

Hi Wolfgang,

You could try something like the routine below, which can be called
after the backwards FFT. This will let each processor store the local
slab of the potential, i.e. the part that's stored locally. It's then
easy to loop over these files in a reading routine and get the full
potential that way.

Volker



void dump_potential(void)
{
  char buf[1000];
  int nprocgroup, masterTask, groupTask, n, i, j, k;
  double asmth, fac, box, tstart, tend;
  float *potential;
  FILE *fd;

  tstart = second();

  if(ThisTask == 0)
    {
      printf("Start dumping potential\n");
      fflush(stdout);
    }

  sprintf(buf, "%s/potential_%d", All.OutputDir, ThisTask);

  nprocgroup = NTask / All.NumFilesWrittenInParallel;

  if((NTask % All.NumFilesWrittenInParallel))
    nprocgroup++;

  masterTask = (ThisTask / nprocgroup) * nprocgroup;

  for(groupTask = 0; groupTask < nprocgroup; groupTask++)
    {
      if(ThisTask == (masterTask + groupTask)) /* ok, it's this
processor's turn */
        {
          if(!(fd = fopen(buf, "w")))
            {
              printf("Error. Can't write in file '%s'\n", buf);
              endrun(11);
            }

          n = PMGRID;
          fwrite(&n, sizeof(int), 1, fd);

          n = sizeof(float);
          fwrite(&n, sizeof(int), 1, fd);

          fwrite(&slabs_per_task[ThisTask], sizeof(int), 1, fd);
          fwrite(&first_slab_of_task[ThisTask], sizeof(int), 1, fd);

          box = All.BoxSize;
          asmth = All.Asmth[0];

          fwrite(&box, sizeof(double), 1, fd);
          fwrite(&asmth, sizeof(double), 1, fd);

          fac = All.G * All.PartMass / (M_PI * All.BoxSize);

          potential = (float *) forcegrid;

          for(i = 0; i < slabs_per_task[ThisTask]; i++)
            for(j = 0; j < PMGRID; j++)
              for(k = 0; k < PMGRID; k++)
                *potential++ = fac * rhogrid[(i * PMGRID + j) * PMGRID2 + k];

          potential = (float *) forcegrid;

          fwrite(potential, sizeof(float), PMGRID * PMGRID *
slabs_per_task[ThisTask], fd);

          fclose(fd);
        }

      /* wait inside the group */
      MPI_Barrier(MPI_COMM_WORLD);
    }


  MPI_Barrier(MPI_COMM_WORLD);

  tend = second();

  if(ThisTask == 0)
    {
      printf("finished writing potential (took=%g sec)\n",
timediff(tstart, tend));
      fflush(stdout);
    }
}
Received on 2005-09-15 16:50:08