ascat.c

Go to the documentation of this file.
00001 /*
00002     This file is part of magpar.
00003 
00004     Copyright (C) 2002-2010 Werner Scholz
00005 
00006     www:   http://www.magpar.net/
00007     email: magpar(at)magpar.net
00008 
00009     magpar is free software; you can redistribute it and/or modify
00010     it under the terms of the GNU General Public License as published by
00011     the Free Software Foundation; either version 2 of the License, or
00012     (at your option) any later version.
00013 
00014     magpar is distributed in the hope that it will be useful,
00015     but WITHOUT ANY WARRANTY; without even the implied warranty of
00016     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
00017     GNU General Public License for more details.
00018 
00019     You should have received a copy of the GNU General Public License
00020     along with magpar; if not, write to the Free Software
00021     Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
00022 */
00023 
00024 static char const Id[] = "$Id: ascat.c 2962 2010-02-04 19:50:44Z scholz $\n\n";
00025 static char const Td[] = "$Today: " __FILE__ "  " __DATE__ "  "  __TIME__  " $\n\n";
00026 
00027 #include "util.h"
00028 
00029 /* Scatters and gathers integer data */
00030 int ascat(int **aindata,int n_items,int s_item,int *newproc,int insize,int send)
00031 {
00032   MagparFunctionInfoBegin;
00033 
00034   int rank,size;
00035   ierr = MPI_Comm_size(PETSC_COMM_WORLD,&size);CHKERRQ(ierr);
00036   ierr = MPI_Comm_rank(PETSC_COMM_WORLD,&rank);CHKERRQ(ierr);
00037 
00038   assert(insize>=1);
00039   /* this function only works for distributing data from first processor
00040      to other processors
00041      TODO: round robin below does not work anyway; thus, simplify this function!
00042   */
00043   if (!rank) {
00044     assert(n_items>0);
00045   }
00046   else {
00047     assert(n_items==0);
00048   }
00049 
00050   int     *ta_n_items;          /* number of items */
00051   ierr = PetscMalloc(insize*sizeof(int),&ta_n_items);CHKERRQ(ierr);
00052 
00053   int     *ta_cnt;              /* counter for items */
00054   ierr = PetscMalloc(insize*sizeof(int),&ta_cnt);CHKERRQ(ierr);
00055 
00056   int     *ta_getn_items;       /* number of items to be received */
00057   ierr = PetscMalloc(insize*sizeof(int),&ta_getn_items);CHKERRQ(ierr);
00058 
00059   int     **ta_procdata;       /* pointer into data arrays */
00060   ierr = PetscMalloc(insize*sizeof(int*),&ta_procdata);CHKERRQ(ierr);
00061 
00062   int     *indata;
00063   indata=*aindata;
00064 
00065   /* allocate memory for temporary array of sorted data;
00066      at least one item to avoid malloc's complaints about allocating zero bytes
00067   */
00068   int *t_data;
00069   ierr = PetscMalloc(n_items*s_item*sizeof(int)+1,&t_data);CHKERRQ(ierr);
00070 
00071   /* reset counters */
00072   for (int i=0; i<insize; i++)
00073   {
00074     ta_n_items[i]=0;
00075     ta_cnt[i]=0;
00076   }
00077 
00078   /* count number of items to be sent */
00079   for (int i=0; i<n_items; i++)
00080   {
00081     assert(newproc[i]>=0 && newproc[i]<insize);
00082     ta_n_items[newproc[i]]++;
00083   }
00084 
00085   /* set pointers to start of data in temporary array */
00086   ta_procdata[0]=t_data;
00087   for (int i=1; i<insize; i++)
00088   {
00089     ta_procdata[i] = ta_procdata[i-1] + s_item*ta_n_items[i-1];
00090   }
00091 
00092   /* check pointers for consistency */
00093   assert(ta_procdata[insize-1]+s_item*ta_n_items[insize-1]==t_data+n_items*s_item);
00094 
00095 
00096   /* copy data sorted into temporary array */
00097   for (int i=0; i<n_items; i++)
00098   {
00099     for (int j=0; j<s_item; j++)
00100     {
00101       *(ta_procdata[newproc[i]]+ta_cnt[newproc[i]]*s_item+j)=indata[i*s_item+j];
00102     }
00103     ta_cnt[newproc[i]]++;
00104   }
00105 
00106   if (!send) {
00107     ierr = PetscFree(*aindata);CHKERRQ(ierr);
00108     *aindata=t_data;
00109 
00110     ierr = PetscFree(ta_procdata);CHKERRQ(ierr);
00111     ierr = PetscFree(ta_n_items);CHKERRQ(ierr);
00112     ierr = PetscFree(ta_cnt);CHKERRQ(ierr);
00113     ierr = PetscFree(ta_getn_items);CHKERRQ(ierr);
00114 
00115     MagparFunctionInfoReturn(0);
00116   }
00117 
00118   /* sendrecv the number of items to be sent/received */
00119   int t_gettot=0;
00120   for (int i=0; i < insize-1; i++) {
00121     int pto,pfrom;
00122     pto   = (rank+i+1+insize) % insize;
00123     pfrom = (rank-i-1+insize) % insize;
00124 
00125 #ifndef UNIPROC
00126     MPI_Status status;
00127     MPI_Sendrecv(&ta_n_items[pto],1,MPI_INT,pto,0,
00128                  &ta_getn_items[pfrom],1,MPI_INT,pfrom,0,PETSC_COMM_WORLD,&status);
00129 #endif
00130 
00131     t_gettot += ta_getn_items[pfrom];
00132   }
00133 
00134   /* allocate array for final local data */
00135   /* at least one byte to avoid failure of malloc: Cannot malloc size zero! */
00136   /* may happen with bndfacvert, that no data are to be delivered */
00137   int *t_fdata;
00138   ierr = PetscMalloc(s_item*(t_gettot+ta_n_items[rank])*sizeof(int)+1,&t_fdata);CHKERRQ(ierr);
00139 
00140   ta_procdata[rank]=t_fdata;
00141 
00142   /* copy local data from original array into final array */
00143   for (int i=0; i<n_items; i++)
00144   {
00145     /* first processor copies its data to beginning of the array */
00146     assert(rank==0);
00147     if (newproc[i]==rank) {
00148       for (int j=0; j<s_item; j++)
00149       {
00150         *ta_procdata[rank]=indata[i*s_item+j];
00151         ta_procdata[rank]++;
00152       }
00153     }
00154   }
00155 
00156   /* check counter for consistency */
00157   assert(ta_procdata[rank]-t_fdata==s_item*ta_n_items[rank]);
00158 
00159   /* sendrecv data to/from all(!) processors */
00160   for (int i=0; i < insize-1; i++) {
00161     int pto,pfrom;
00162     pto   = (rank+i+1+insize) % insize;
00163     pfrom = (rank-i-1+insize) % insize;
00164 
00165 #ifndef UNIPROC
00166     MPI_Status status;
00167     MPI_Sendrecv(ta_procdata[pto],s_item*ta_n_items[pto],MPI_INT,pto,1,
00168                  ta_procdata[rank],s_item*ta_getn_items[pfrom],MPI_INT,pfrom,1,PETSC_COMM_WORLD,&status);
00169 #endif
00170 
00171     ta_procdata[rank] += s_item*ta_getn_items[pfrom];
00172   }
00173 
00174   /* check counter for consistency */
00175   assert(ta_procdata[rank]-t_fdata==s_item*(t_gettot+ta_n_items[rank]));
00176 
00177   ierr = PetscFree(t_data);CHKERRQ(ierr);
00178   if (n_items>0) {
00179     ierr = PetscFree(*aindata);CHKERRQ(ierr);
00180   }
00181   *aindata=t_fdata;
00182 
00183   t_gettot += ta_n_items[rank];
00184 
00185   ierr = PetscFree(ta_procdata);CHKERRQ(ierr);
00186   ierr = PetscFree(ta_n_items);CHKERRQ(ierr);
00187   ierr = PetscFree(ta_cnt);CHKERRQ(ierr);
00188   ierr = PetscFree(ta_getn_items);CHKERRQ(ierr);
00189 
00190   MagparFunctionInfoReturn(t_gettot);
00191 }
00192 

magpar - Parallel Finite Element Micromagnetics Package
Copyright (C) 2002-2010 Werner Scholz