path: root/src/SetupPGH.c

 /*@@
   @file      SetupPGH.c
   @date      Fri Feb 21 10:18:17 1997
   @author    Paul Walker
   @desc 
   Initializes the Pugh Grid Hierachy. 
   @enddesc 
   @version $Header$
 @@*/

/*#define DEBUG_PUGH*/

#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <assert.h>

#ifdef HAVE_UNISTD_H
#include <unistd.h>
#endif

#include "cctk.h"
#include "pugh.h"
#include "pughi.h"
#include "cctk_Parameters.h"

static char *rcsid = "$Header$";

CCTK_FILEVERSION(CactusPUGH_PUGH_SetupPGH_c)



 /*@@
   @routine    PUGH_SetupPGH
   @date       Fri Feb 21 10:21:36 1997
   @author     Gabrielle Allen, Tom Goodale, Paul Walker
   @desc 
   Sets up the PGH distribution based on processor number and
   problem size. This includes setting an optimal domain
   decomposition for each dimension, and setting ghosts and overlaps.
   @enddesc 
   @hdate 4 June 1999
   @hauthor Thomas Radke
   @hdesc Check for unique precision of CCTK_INT and CCTK_REAL types
          in metacomputing environment
   @hdate 30 June 1999
@@*/

pGH *PUGH_SetupPGH(void *callerid, 
                   int dim, 
                   int *nsize, 
                   int *nghostzones,
                   int staggertype,
                   int *perme) 
{
  DECLARE_CCTK_PARAMETERS

  pGH *GH;
  int idim;


  /* Allocate for myself */
  GH = (pGH*) malloc(sizeof(pGH));

  /* Set an identifier for my parent */
  GH->callerid = callerid ;

  GH->GFExtras     = (pGExtras **)malloc(dim*sizeof(pGExtras*));
  GH->Connectivity = (pConnectivity **)malloc(dim*sizeof(pConnectivity*));

  /* Set the total number of processors */
  pGH_SetupnProcs(GH,dim);

  /* Set up connectivity and extras for each dimension */
  for (idim=1;idim<=dim;idim++)
  {
    int *nprocs;

    nprocs = (int *) malloc(idim*sizeof(int));

    PUGH_SetupDefaultTopology(idim,nprocs);

    GH->Connectivity[idim-1]=
      PUGH_SetupConnectivity(idim, 
                             GH->nprocs, 
                             nprocs,
                             perme);
    free(nprocs);
      
    GH->GFExtras[idim-1]=
      PUGH_SetupPGExtras(idim, 
                         perme,
                         staggertype,
                         nsize,
                         nghostzones,
                         GH->nprocs,
                         ((pConnectivity **)GH->Connectivity)[idim-1]->nprocs,
                         GH->myproc);

  }

  GH->active = 1;
  GH->commmodel = PUGH_ALLOCATEDBUFFERS;

  /* create the timer for communication time */
  if(timer_output)
  {
    GH->comm_time = CCTK_TimerCreateI ();
  }
  else
  {
    GH->comm_time = -1;
  }

  /* set staggering flag */
  /*GH->stagger = staggertype;                 */
  GH->identity_string = NULL;
  GH->level = 0;
  GH->mglevel = 0;
  GH->convlevel = 0;
  GH->forceSync = 0;
  GH->nvariables = 0;
  GH->narrays = 0;
  GH->variables = NULL;
  return GH;

}


int pGH_SetupnProcs(pGH *GH, int dim)
{
#ifdef CCTK_MPI
  CCTK_REAL4 sum_sizes [2], compiled_sizes [2];
#endif 

  DECLARE_CCTK_PARAMETERS;

  GH->dim = dim;

#ifdef CCTK_MPI
  /* Set up my communicator. This would allow us to run pugh on 
     a subset of processors at a later date if, for instance, we
     were using panda or what not.
   */
  CACTUS_MPI_ERROR(MPI_Comm_dup(MPI_COMM_WORLD, &(GH->PUGH_COMM_WORLD)));

  CACTUS_MPI_ERROR(MPI_Comm_size(GH->PUGH_COMM_WORLD, &(GH->nprocs)));
  CACTUS_MPI_ERROR(MPI_Comm_rank(GH->PUGH_COMM_WORLD, &(GH->myproc)));

  /* check that all executables uses the same integer and fp precision
     within a metacomputing environment
     NOTE: We cannot use CCTK_INT4s in MPI_Allreduce() since
           (although they have the same size) they might be correspond to
           different MPI datatypes.
           CCTK_REAL4 should always refer to MPI_FLOAT. */
  compiled_sizes [0] = sizeof (CCTK_INT);
  compiled_sizes [1] = sizeof (CCTK_REAL);
  CACTUS_MPI_ERROR (MPI_Allreduce (compiled_sizes, sum_sizes, 2, PUGH_MPI_REAL4,
                    MPI_SUM, GH->PUGH_COMM_WORLD));
  if (compiled_sizes [0] * GH->nprocs != sum_sizes [0])
  {
    CCTK_WARN (0, "Cannot run executables with different precision for "
                  "CCTK_INTs within a metacomputing environment !\n"
                  "Please configure with unique CCTK_INTEGER_PRECISION !");
  }

  if (compiled_sizes [1] * GH->nprocs != sum_sizes [1])
  {
    CCTK_WARN (0, "Cannot run executables with different precision for "
                  "CCTK_REALs within a metacomputing environment !\n"
                  "Please configure with unique CCTK_REAL_PRECISION !");
  }

  /* define the complex datatype as a concatanation of 2 PUGH_MPI_REALs */
  CACTUS_MPI_ERROR(MPI_Type_contiguous (2, PUGH_MPI_REAL, &GH->PUGH_mpi_complex));
  CACTUS_MPI_ERROR(MPI_Type_commit (&GH->PUGH_mpi_complex));
#else
  GH->nprocs = 1;
  GH->myproc = 0;
#endif

  return 0;
}


int PUGH_Terminate (cGH *GH)
{
  pGH *pughGH = PUGH_pGH (GH);


  PUGH_DestroyPGH (&pughGH);
  return (0);
}


 /*@@
   @routine    PUGH_DestroyPGH
   @date       Thu Aug 21 11:38:10 1997
   @author     Paul Walker
   @desc 
   Destroys a GH object.
   @enddesc 
@@*/

void PUGH_DestroyPGH(pGH **GHin) 
{
  /* First remove me from the list. */
  pGH    *GH;
  pGA    *GA;
  cGroup pgroup;
  int i;
  int variable;
  int group;
  int this_var;

  GH = *GHin;

#ifdef CCTK_MPI
  CACTUS_MPI_ERROR(MPI_Type_free (&GH->PUGH_mpi_complex));
  CACTUS_MPI_ERROR(MPI_Comm_free(&(GH->PUGH_COMM_WORLD)));
#endif

  /* Great. Now go about the work of destroying me. */

  variable = 0;

  for(group = 0; group < CCTK_NumGroups(); group++)
  {

#ifdef DEBUG_PUGH
    printf("Calling Destroying Group %s\n", CCTK_GroupName(group));
    fflush(stdout);
#endif 

    CCTK_GroupData(group,&pgroup);

    /* destroy group comm buffers */
    if (pgroup.grouptype == CCTK_ARRAY || pgroup.grouptype == CCTK_GF)
    {
      GA = (pGA *) GH->variables[variable][0];

      /* Destroy group comm buffers */
      if (GA->groupcomm)
      {
        if (GA->groupcomm->commflag != PUGH_NOCOMM)
        {
          PUGH_DisableGArrayGroupComm (GH, variable, GA->groupcomm);
        }
        PUGH_DestroyComm (&GA->groupcomm);
      }

      /* Destroy the group's connectivity and extras structure
         for CCTK_ARRAY groups.
         Remember that the connectivity and extras for CCTK_GF types
         are shared between all such groups and are destroyed later. */
      if (GA->connectivity != GH->Connectivity[pgroup.dim-1])
      {
        PUGH_DestroyConnectivity (&GA->connectivity);
      }
      if (GA->extras != GH->GFExtras[pgroup.dim-1])
      {
        PUGH_DestroyPGExtras (&GA->extras);
      }
    }

    for (this_var = 0; this_var < pgroup.numvars; this_var++, variable++) 
    {
      for(i = 0 ; i < pgroup.numtimelevels; i++)
      {
        switch(pgroup.grouptype)
        {
          case CCTK_GF:
          case CCTK_ARRAY:
            PUGH_DestroyGArray(&(((pGA ***)GH->variables)[variable][i]));
            break;
          case CCTK_SCALAR:
            free(GH->variables[variable][i]);
            break;
        }
      }
      free(GH->variables[variable]);
    }
  }

  for (i=1;i<=GH->dim;i++)
  {
    PUGH_DestroyConnectivity(&GH->Connectivity[i-1]);
    PUGH_DestroyPGExtras(&GH->GFExtras[i-1]);
  }

  if(GH->identity_string)
  {
    free(GH->identity_string);
  }
  free(GH->Connectivity);
  free(GH->GFExtras);
  if(GH->variables)
  {
    free(GH->variables);
  }
  free(GH);
  *GHin = NULL;
}

void pGH_DumpInfo(pGH *GH) 
{
  int i,j,k;
  printf("INFO: \n");
  printf("   myproc: %d/%d \n",GH->myproc,GH->nprocs);
  printf("   ownership: \n");
  for (i=0;i<PUGH_NSTAGGER;i++)
  {
    for (j=0;j<GH->dim;j++)
    {
      printf("     GH->ownership[%d][0/1][%d]: %d %d \n",
             i,j,GH->GFExtras[GH->dim-1]->ownership[i][0][j],GH->GFExtras[GH->dim-1]->ownership[i][1][j]);
    }
  }

  for (i=0;i<PUGH_NSTAGGER;i++)
  {
    for (j=0;j<GH->dim;j++)
    {
      for (k=0;k<2*GH->dim;k++)
      {
        printf("     GH->GFExtras[GH->dim-1]->ghosts[%d][0/1][%d][%d]: %d %d \n",
               i,k,j,GH->GFExtras[GH->dim-1]->ghosts[i][0][k][j],GH->GFExtras[GH->dim-1]->ghosts[i][1][k][j]);
      }
    }
  }

  for (i=0;i<PUGH_NSTAGGER;i++)
  {
    for (j=0;j<GH->dim;j++)
    {
      for (k=0;k<2*GH->dim;k++)
      {
        printf("     GH->GFExtras[GH->dim-1]->overlap[%d][0/1][%d][%d]: %d %d \n",
               i,k,j,GH->GFExtras[GH->dim-1]->overlap[i][0][j][k],GH->GFExtras[GH->dim-1]->overlap[i][1][k][j]);
      }
    }
  }
}

 /*@@
   @routine    PUGH_SetupDefaultTopology
   @date       Sun 23 Jan
   @author     Gabrielle Allen
   @desc 
   Use PUGH parameters to set the processor decomposition 
   is required
   @enddesc 
   @hdate 
   @hauthor 
   @hdesc 
   @hdate 
@@*/

int PUGH_SetupDefaultTopology(int dim, int *nprocs)
{

  DECLARE_CCTK_PARAMETERS

  int retval=PUGH_SUCCESS;

  switch(dim)
  {
    case 1:
      nprocs[0] = processor_topology_1d_x;
      break;
    case 2:
      nprocs[0] = processor_topology_2d_x;
      nprocs[1] = processor_topology_2d_y;
      break;
    case 3:
      nprocs[0] = processor_topology_3d_x;
      nprocs[1] = processor_topology_3d_y;
      nprocs[2] = processor_topology_3d_z;
      break;
    default:
      nprocs=NULL;
      retval=PUGH_ERROR;
  }

  return retval;
}


int PUGH_ParallelInit(cGH *GH)
{
  return 0;
}

int PUGH_Exit(cGH *GH, int retval)
{
#ifdef CCTK_MPI
  CACTUS_MPI_ERROR(MPI_Finalize());
#endif
  exit(retval);
  return(retval);
}

int PUGH_Abort(cGH *GH, int retval)
{
#ifdef CCTK_MPI
  MPI_Comm COMM;
#endif

#ifdef CCTK_MPI
  if (GH)
  {
    COMM = PUGH_pGH(GH)->PUGH_COMM_WORLD;
  }
  else
  {
    COMM = MPI_COMM_WORLD;
  }
  CACTUS_MPI_ERROR(MPI_Abort(COMM,retval));
#else
  /* FIXME */
  /*abort();*/
#endif
  exit(0);
  return(0);
}


int PUGH_MyProc(cGH *GH)
{
  int myproc;

#ifdef CCTK_MPI
  MPI_Comm comm;
  if(GH)
  {
    comm = PUGH_pGH(GH)->PUGH_COMM_WORLD;
  }
  else
  {
    comm = MPI_COMM_WORLD;
  }
  CACTUS_MPI_ERROR(MPI_Comm_rank(comm, &myproc));
#else
  myproc = 0;
#endif

  return myproc;
}


int PUGH_nProcs(cGH *GH)
{
  int nprocs;

#ifdef CCTK_MPI
  MPI_Comm comm;
  if(GH)
  {
    comm = PUGH_pGH(GH)->PUGH_COMM_WORLD;
  }
  else
  {
    comm = MPI_COMM_WORLD;
  }
  CACTUS_MPI_ERROR(MPI_Comm_size(comm, &nprocs));
#else
  nprocs = 1;
#endif

  return nprocs;
}