src/kernel/calcgrid.c

   1 /*
   2  * $Id$
   3  *
   4  *       This source code is part of
   5  *
   6  *        G   R   O   M   A   C   S
   7  *
   8  * GROningen MAchine for Chemical Simulations
   9  *
  10  *               VERSION 2.0
  11  *
  12  * Copyright (c) 1991-1999
  13  * BIOSON Research Institute, Dept. of Biophysical Chemistry
  14  * University of Groningen, The Netherlands
  15  *
  16  * Please refer to:
  17  * GROMACS: A message-passing parallel molecular dynamics implementation
  18  * H.J.C. Berendsen, D. van der Spoel and R. van Drunen
  19  * Comp. Phys. Comm. 91, 43-56 (1995)
  20  *
  21  * Also check out our WWW page:
  22  * http://md.chem.rug.nl/~gmx
  23  * or e-mail to:
  24  * gromacs@chem.rug.nl
  25  *
  26  * And Hey:
  27  * GRowing Old MAkes el Chrono Sweat
  28  */
  29 static char *SRCID_calcgrid_c = "$Id$";
  30
  31 #include "typedefs.h"
  32 #include "smalloc.h"
  33 #include "fatal.h"
  34 #include "calcgrid.h"
  35
  36 #define facNR 6
  37 int factor[facNR] = {2,3,5,7,11,13};
  38 int decomp[facNR];
  39 int ng,ng_max,*list,n_list,n_list_alloc;
  40
  41 static void make_list(int start_fac)
  42 {
  43   int i;
  44
  45   if (ng < ng_max) {
  46     if (n_list >= n_list_alloc) {
  47       n_list_alloc += 100;
  48       srenew(list,n_list_alloc);
  49     }
  50     list[n_list] = ng;
  51     n_list++;
  52
  53     for(i=start_fac; i<facNR; i++) {
  54       /* allow any power of 2, 3, 5 and 7, but only one of 11 or 13 */
  55       if (i<4 || (decomp[4]+decomp[5]==0)) {
  56         ng*=factor[i];
  57         decomp[i]++;
  58         make_list(i);
  59         ng/=factor[i];
  60         decomp[i]--;
  61       }
  62     }
  63   }
  64 }
  65
  66 static int list_comp(const void *a,const void *b)
  67 {
  68   return (*((int *)a) - *((int *)b));
  69 }
  70
  71 real calc_grid(matrix box,real gr_sp,int *nx,int *ny,int *nz,int nprocs)
  72 {
  73   int  d,n[DIM];
  74   int  i,nmin[DIM];
  75   rvec box_size,spacing;
  76   real max_spacing;
  77
  78   if (gr_sp <= 0)
  79     fatal_error(0,"invalid fourier grid spacing: %g",gr_sp);
  80
  81   for(d=0; d<DIM; d++)
  82     box_size[d] = box[d][d];
  83
  84   n[XX] = *nx;
  85   n[YY] = *ny;
  86   n[ZZ] = *nz;
  87
  88   ng = 1;
  89   ng_max = 1;
  90   for(d=0; d<DIM; d++) {
  91     nmin[d] = (int)(box_size[d]/gr_sp + 0.999);
  92     if (2*nmin[d] > ng_max)
  93       ng_max = 2*nmin[d];
  94   }
  95   n_list=0;
  96   n_list_alloc=0;
  97   list=NULL;
  98   for(i=0; i<facNR; i++)
  99     decomp[i]=0;
 100   make_list(0);
 101
 102   if ((*nx<=0) || (*ny<=0) || (*nz<=0))
 103     fprintf(stderr,"Calculating fourier grid dimensions for%s%s%s\n",
 104             *nx > 0 ? "":" X",*ny > 0 ? "":" Y",*nz > 0 ? "":" Z");
 105
 106   qsort(list,n_list,sizeof(list[0]),list_comp);
 107   if (debug)
 108     for(i=0; i<n_list; i++)
 109       fprintf(debug,"grid: %d\n",list[i]);
 110
 111   if (((*nx>0) && (*nx != nprocs*(*nx/nprocs))) ||
 112       ((*ny>0) && (*ny != nprocs*(*ny/nprocs))))
 113     fatal_error(0,"the x or y grid spacing (nx %d, ny %d) is not divisible by the number of processors (%d)",*nx,*ny,nprocs);
 114
 115   for(d=0; d<DIM; d++) {
 116     for(i=0; (i<n_list) && (n[d]<=0); i++)
 117       if ((list[i] >= nmin[d]) &&
 118           ((d == ZZ) || (list[i] == nprocs*(list[i]/nprocs))))
 119         n[d] = list[i];
 120     if (n[d] <= 0)
 121       fatal_error(0 ,"could not find a grid spacing with nx and ny divisible by the number of processors (%d)",nprocs);
 122   }
 123
 124   max_spacing = 0;
 125   for(d=0; d<DIM; d++) {
 126     spacing[d] = box_size[d]/n[d];
 127     if (spacing[d] > max_spacing)
 128       max_spacing = spacing[d];
 129   }
 130   *nx = n[XX];
 131   *ny = n[YY];
 132   *nz = n[ZZ];
 133   fprintf(stderr,"Using a fourier grid of %dx%dx%d, spacing %.3f %.3f %.3f\n",
 134           *nx,*ny,*nz,spacing[XX],spacing[YY],spacing[ZZ]);
 135
 136   return max_spacing;
 137 }
 138
 139