diff options
| -rw-r--r-- | README | 7 | ||||
| -rw-r--r-- | interface.ccl | 5 | ||||
| -rw-r--r-- | param.ccl | 2 | ||||
| -rw-r--r-- | schedule.ccl | 8 | ||||
| -rw-r--r-- | src/CollectData1D.c | 279 | ||||
| -rw-r--r-- | src/CollectData2D.c | 47 | ||||
| -rw-r--r-- | src/GetLocalLine.c | 336 | ||||
| -rw-r--r-- | src/TestSlab.c | 103 | ||||
| -rw-r--r-- | src/make.code.defn | 9 |
9 files changed, 796 insertions, 0 deletions
@@ -0,0 +1,7 @@ +Cactus Code Thorn Hyperslab +Authors : ... +CVS info : $Header$ +-------------------------------------------------------------------------- + +Purpose of the thorn: + diff --git a/interface.ccl b/interface.ccl new file mode 100644 index 0000000..9f47f5b --- /dev/null +++ b/interface.ccl @@ -0,0 +1,5 @@ +# Interface definition for thorn Hyperslab +# $Header$ + + +implements: hyperslab diff --git a/param.ccl b/param.ccl new file mode 100644 index 0000000..1a19316 --- /dev/null +++ b/param.ccl @@ -0,0 +1,2 @@ +# Parameter definitions for thorn Hyperslab +# $Header$ diff --git a/schedule.ccl b/schedule.ccl new file mode 100644 index 0000000..e35b24c --- /dev/null +++ b/schedule.ccl @@ -0,0 +1,8 @@ +# Schedule definitions for thorn Hyperslab +# $Header$ + +schedule TestSlab at INITIAL +{ + LANG:C +} "Testing slabbing" + diff --git a/src/CollectData1D.c b/src/CollectData1D.c new file mode 100644 index 0000000..c35fb59 --- /dev/null +++ b/src/CollectData1D.c @@ -0,0 +1,279 @@ + +#include <stdio.h> +#include <stdlib.h> +#include <math.h> + +#include "cctk.h" +#include "cctk_Parameters.h" +#include "cctk_Arguments.h" + +#include "CactusBase/IOUtil/src/ioGH.h" +#include "CactusPUGH/PUGH/src/include/pugh.h" + +/* MPI message tag */ +#define TAGBASE 1234 +#define MAX_DIM 3 + +int CollectData1D(cGH *GH, int vindex, int vtype, int vtypesize,int vtimelvl, + int vdim, int *S_l, int *u_l, int ds, int len_l, + int dsize, void **dptr) +{ + + void *loc_dptr; + int loc_dcount=0; + + void *all_dptr; + int all_dcount=0; + + int *rem_dcount; + + int iproc,myproc,nprocs; + int vmpi_type; + int *displs; + + int idim, i,j,k, ip; + int ierr; + + pGH *pughGH; + + pughGH = PUGH_pGH (GH); + myproc = CCTK_MyProc(GH); + nprocs = CCTK_nProcs(GH); + vmpi_type= PUGH_MPI_REAL; + + rem_dcount=(CCTK_INT *) malloc (nprocs * sizeof (CCTK_INT)); + + ierr = CollectLocalData1D(GH, vindex, vtype, vtypesize, vtimelvl, + vdim, S_l, u_l, len_l, + &loc_dptr, &loc_dcount); + + + /* Send local number of data elements */ + CACTUS_MPI_ERROR (MPI_Gather (&loc_dcount, 1, vmpi_type, + rem_dcount, 1, PUGH_MPI_INT, 0, pughGH->PUGH_COMM_WORLD)); + + + /* Calculate the displacements */ + displs[0] = 0; + all_dcount= rem_dcount[0]; + for (iproc=1;iproc<nprocs;iproc++) + { + displs[iproc] = displs[iproc-1] + rem_dcount[iproc-1]; + all_dcount += rem_dcount[ip]; + } + + all_dptr = malloc (vtypesize*all_dcount); + + CACTUS_MPI_ERROR (MPI_Gatherv(loc_dptr, loc_dcount, vmpi_type, + all_dptr, rem_dcount, displs, vmpi_type, + 0, pughGH->PUGH_COMM_WORLD)); + + return(0); + + +} + + +int CollectLocalData1D(cGH *GH, int vindex, int vtype, int vtypesize, int vtimelvl, + int vdim, int *S_l, int *u_l, int ds, int len_l, + void **loc_dptr, int *dsize) +{ + int *locstart, *locend, nlocpoints; + int ipnt, gridpnt[MAX_DIM]; + int ierr; + + int idim,l; + int vmpi_type; + + void *data = CCTK_VarDataPtrI (GH, vtimelvl, vindex); + + ioGH *ioUtilGH; + + ioUtilGH = (ioGH *) GH->extensions [CCTK_GHExtensionHandle ("IO")]; + + + switch (vtype) { + + case CCTK_VARIABLE_CHAR: + vtypesize = sizeof (CCTK_CHAR); + vmpi_type = PUGH_MPI_CHAR; + break; + + case CCTK_VARIABLE_INT: + vtypesize = sizeof (CCTK_INT); + vmpi_type = PUGH_MPI_INT; + break; + + case CCTK_VARIABLE_REAL: + vtypesize = ioUtilGH->out_single ? sizeof(CCTK_REAL4):sizeof (CCTK_REAL); + vmpi_type = ioUtilGH->out_single ? PUGH_MPI_REAL4 : PUGH_MPI_REAL; + break; + + default: + CCTK_WARN (1, "Unsupported variable type in IOFlexIO_Write2D"); + return(-1); + } + + locstart = (int*)malloc(vdim*sizeof(int)); + locend = (int*)malloc(vdim*sizeof(int)); + + ierr = GetLocalLine(GH, vdim, S_l, u_l, ds, len_l, locstart, locend, &nlocpoints); + + printf("CollectLocal1D: %d npoints\n",nlocpoints); + + /*$ *loc_dptr= malloc(nlocpoints * vtypesize);$*/ + + l=0; + + /* CCTK_VARIABLE_INT */ + if (vtype == CCTK_VARIABLE_INT) + { + CCTK_INT *int_dptr = (CCTK_INT *) *loc_dptr; + + for (ipnt=0; ipnt<=nlocpoints; ipnt++) + { + + if (vdim==3) { + for (ipnt=0; ipnt<=nlocpoints; ipnt++) + { + for (idim=0;idim<vdim;idim++) + gridpnt[idim] = locstart[idim]+ipnt * u_l[idim]; + int_dptr[l++]=((CCTK_INT *) data) + [CCTK_GFINDEX3D(GH,gridpnt[0],gridpnt[1],gridpnt[2])]; + } + } + else if (vdim==2) { + for (ipnt=0; ipnt<=nlocpoints; ipnt++) + { + for (idim=0;idim<vdim;idim++) + gridpnt[idim] = locstart[idim]+ipnt * u_l[idim]; + int_dptr[l++]=((CCTK_INT *) data) + [CCTK_GFINDEX2D(GH,gridpnt[0],gridpnt[1])]; + } + } + else if (vdim==1) { + for (ipnt=0; ipnt<=nlocpoints; ipnt++) + { + for (idim=0;idim<vdim;idim++) + gridpnt[idim] = locstart[idim]+ipnt * u_l[idim]; + int_dptr[l++]=((CCTK_INT *) data) + [CCTK_GFINDEX1D(GH,gridpnt[0])]; + } + } + else return(-1); + } + } + /* CCTK_VARIABLE_CHAR */ + else if (vtype == CCTK_VARIABLE_CHAR) + { + CCTK_CHAR *char_dptr = (CCTK_CHAR *) *loc_dptr; + + if (vdim==3) { + for (ipnt=0; ipnt<=nlocpoints; ipnt++) + { + for (idim=0;idim<vdim;idim++) + gridpnt[idim] = locstart[idim]+ipnt * u_l[idim]; + char_dptr[l++]=((CCTK_CHAR *) data)[CCTK_GFINDEX3D(GH,gridpnt[0],gridpnt[1],gridpnt[2])]; + } + } + else if (vdim==2) { + for (ipnt=0; ipnt<=nlocpoints; ipnt++) + { + for (idim=0;idim<vdim;idim++) + gridpnt[idim] = locstart[idim]+ipnt * u_l[idim]; + char_dptr[l++]=((CCTK_CHAR *) data) + [CCTK_GFINDEX2D(GH,gridpnt[0],gridpnt[1])]; + } + } + else if (vdim==1) { + for (ipnt=0; ipnt<=nlocpoints; ipnt++) + { + for (idim=0;idim<vdim;idim++) + gridpnt[idim] = locstart[idim]+ipnt * u_l[idim]; + char_dptr[l++]=((CCTK_CHAR *) data) + [CCTK_GFINDEX1D(GH,gridpnt[0])]; + } + } + } + + /* CCTK_VARIABLE_REAL */ + else if (vtype == CCTK_VARIABLE_REAL && !ioUtilGH->out_single) + { + CCTK_REAL *real_dptr = (CCTK_REAL *) *loc_dptr; + + + if (vdim==3) { + for (ipnt=0; ipnt<=nlocpoints; ipnt++) + { + for (idim=0;idim<vdim;idim++) + gridpnt[idim] = locstart[idim]+ipnt * u_l[idim]; + real_dptr[l++]=((CCTK_REAL *) data) + [CCTK_GFINDEX3D(GH,gridpnt[0],gridpnt[1],gridpnt[2])]; + } + } + else if (vdim==2) { + for (ipnt=0; ipnt<=nlocpoints; ipnt++) + { + for (idim=0;idim<vdim;idim++) + gridpnt[idim] = locstart[idim]+ipnt * u_l[idim]; + real_dptr[l++]=((CCTK_REAL *) data) + [CCTK_GFINDEX2D(GH,gridpnt[0],gridpnt[1])]; + } + } + else if (vdim==1) { + for (ipnt=0; ipnt<=nlocpoints; ipnt++) + { + for (idim=0;idim<vdim;idim++) + gridpnt[idim] = locstart[idim]+ipnt * u_l[idim]; + real_dptr[l++]=((CCTK_REAL *) data) + [CCTK_GFINDEX1D(GH,gridpnt[0])]; + } + } + } + /* CCTK_VARIABLE_REAL4 */ + else if (vtype == CCTK_VARIABLE_REAL4 && ioUtilGH->out_single) + { + CCTK_REAL4 *real4_dptr = (CCTK_REAL4*) *loc_dptr; + + + if (vdim==3) { + for (ipnt=0; ipnt<=nlocpoints; ipnt++) + { + for (idim=0;idim<vdim;idim++) + gridpnt[idim] = locstart[idim]+ipnt * u_l[idim]; + + real4_dptr[l++]=((CCTK_REAL4 *) data) + [CCTK_GFINDEX3D(GH,gridpnt[0],gridpnt[1],gridpnt[2])]; + } + } + else if (vdim==2) { + for (ipnt=0; ipnt<=nlocpoints; ipnt++) + { + for (idim=0;idim<vdim;idim++) + gridpnt[idim] = locstart[idim]+ipnt * u_l[idim]; + + real4_dptr[l++]=((CCTK_REAL4 *) data) + [CCTK_GFINDEX2D(GH,gridpnt[0],gridpnt[1])]; + } + } + else if (vdim==1) { + for (ipnt=0; ipnt<=nlocpoints; ipnt++) + { + for (idim=0;idim<vdim;idim++) + gridpnt[idim] = locstart[idim]+ipnt * u_l[idim]; + + real4_dptr[l++]=((CCTK_REAL4 *) data) + [CCTK_GFINDEX1D(GH,gridpnt[0])]; + } + } + } +} + + + + + + + + + diff --git a/src/CollectData2D.c b/src/CollectData2D.c new file mode 100644 index 0000000..4da603c --- /dev/null +++ b/src/CollectData2D.c @@ -0,0 +1,47 @@ + +#include <stdio.h> +#include <stdlib.h> +#include <math.h> + +#include "cctk.h" +#include "cctk_Parameters.h" +#include "cctk_Arguments.h" + +/* MPI message tag */ +#define TAGBASE 1234 + +int CollectData2D(cGH *GH, int vindex, int vtype, int vtypesize,int vtimelvl, + int vdim, int *S_l, int *u1_l, int *u2_l, int *len_l, + int *dsize, void **dptr) +{ + int ierr; + + ierr = CollectLocalData2D(GH, vindex, vtype, vtypesize, vtimelvl, + vdim, S_l, u1_l, u2_l, len_l, + &loc_dptr, &loc_dcount); + +} + +int CollectLocalData2D(cGH *GH, int vindex, int vtype, int vtypesize, int vtimelvl, + int vdim, int *S_l, int *u1_l,int *u2_l, int *len_l, + void **loc_dptr, int *dsize) + +{ + + int idim; + int t1,t2,t3; + int gridpnt[3]; + + + ierr = GetLocalLine(GH, vdim, S_l, u1_l, len_l[0], locstart1, locend1, &nlocpoints1); + ierr = GetLocalLine(GH, vdim, S_l, u1_l, len_l[0], locstart2, locend2, &nlocpoints2); + + for (l1=0;l1<nlocpoints1;l1++) + for (l2=0;l2<nlocpoints2;l2++) + { + for (idim=0;idim<vdim;idim++) + gridpnt[idim]=S_l[idim]+l1*u1_l[idim]+l2*u2_l[idim]; + real_ptr[l++] = ((CCTK_REAL*) data)[CCTK_GFINDEX3D(GH,gridpnt[0],gridpnt[1],gridpnt[2])]; + } + +} diff --git a/src/GetLocalLine.c b/src/GetLocalLine.c new file mode 100644 index 0000000..c0b1ef3 --- /dev/null +++ b/src/GetLocalLine.c @@ -0,0 +1,336 @@ + + +/* + A B + <-----------------> + |_____|_____|_____|_____|_____|_____|_____| + 0 1 2 3 4 5 6 #proc + + These are the case we have to deal with. For downsampling we need + to calulate the offset that a proc needs to start on gridpoint, that + is actually reached by downsampling. Eg. see the "X" below. + where downsampl. = 4, start at global index 2, 3 procs, shown with + ghostzones, where they overlap. + 1 1 1 1 1 1 1 1 1 1 + globalIDX: 0 1 3 4 5 6 7 8 9 0 1 2 3 4 5 6 7 8 9 + + localIDX: + proc0/2: 0 1 2 3 4 5 6 7 0 1 2 3 4 5 6 + proc1: 0 1 2 3 4 5 6 7 + downsampl: X X X X + All the stuff below is about starting at "4" on proc one and + not 0 or 1 or what else. (0 on proc 2 is a ghost). + +*/ + + +#include <stdio.h> +#include <stdlib.h> +#include <math.h> + +#include "cctk.h" +#include "cctk_Parameters.h" +#include "cctk_Arguments.h" + +#define MAX_DIM 3 +#define MAX_FACE 6 +#define EPS 1e-10 +#define BAD -42 + +#define ABS(a) ((a)<0 ? -(a) : (a)) +#define SIGN(a) ((a)<0.0 ? (-1.0) : (1.0)) + +int IsPointOnProc_INT (int dim, int *gownership, int *pnt); +int IsPointOnProc_REAL(int dim, int *gownership, CCTK_REAL *pnt); + +int IsPointOnProc_INT(int dim, int *gownership, int *pnt) +{ + int retval=1; /*assume YES */ + int idim; + for (idim=0;idim<dim;idim++) + { + /* out of range -> set retval to zero;*/ + if ((gownership[2*idim ]>pnt[idim]) || + (gownership[2*idim+1]<pnt[idim])) + retval&=0; /* incremental AND on the return value */ + } + return(retval); +} + +int IsPointOnProc_REAL(int dim, int *gownership, CCTK_REAL *pnt) +{ + int retval=1; /*assume YES */ + int idim; + for (idim=0;idim<dim;idim++) + { + /* out of range -> set retval to zero;*/ + if ((gownership[2*idim ]>pnt[idim]) || + (gownership[2*idim+1]<pnt[idim])) + retval&=0; /* incremental AND on the return value */ + } + return(retval); +} + +int GetLocalLine(cGH *GH, + int vardim, int *S_l, int *u_l, int *ds, int len_l, + int *locstart, int *locend, int *nlocpoints) +{ + /* variables for proc. surfaces */ + CCTK_REAL lmd_s[MAX_FACE], lmd[2], t_lmd; + int P_s[MAX_DIM][MAX_FACE]; + int n_s[MAX_DIM][MAX_FACE]; + int E_l[MAX_DIM]; + int notcoll[MAX_FACE],isloc[2]; + CCTK_REAL t_pnt1[MAX_DIM], t_pnt2[MAX_DIM]; + CCTK_REAL secpnt[MAX_DIM][2]; + int isec; + + int prochasline=0; + + int tmp1,tmp2,t; + CCTK_REAL rtmp1,rtmp2; + + int idim, iface, myproc; + int *gownership; + + myproc = CCTK_MyProc(GH); + + gownership= (int*) malloc(2*vardim*sizeof(int)); + + printf("\n+++++++++\n"); + printf("S_l [XYZ] : %d %d %d \n", + S_l[0], S_l[1], S_l[2]); + printf("u_l [XYZ] : %d %d %d \n", + u_l[0], u_l[1], u_l[2]); + + for (idim=0;idim<vardim;idim++) + { + /* GLOBAL OWNERSHIP: lower/upperbnd are INCLUSIVE */ + gownership[2*idim] = GH->cctk_lbnd[idim]; + if (GH->cctk_bbox[2*idim]==0) + gownership[2*idim] += GH->cctk_nghostzones[idim]; + + gownership[2*idim+1]= GH->cctk_ubnd[idim]; + if (GH->cctk_bbox[2*idim+1]==0) + gownership[2*idim+1]-= (GH->cctk_nghostzones[idim]); + } + + + for (idim=0;idim<vardim;idim++) + { + P_s[idim][0] = gownership[2*idim]; /* Aufpunkt, minus faces */ + P_s[idim][2] = gownership[2*idim]; /* Aufpunkt, minus faces */ + P_s[idim][4] = gownership[2*idim]; /* Aufpunkt, minus faces */ + P_s[idim][1] = gownership[2*idim+1]; /* Aufpunkt, plus faces */ + P_s[idim][3] = gownership[2*idim+1]; /* Aufpunkt, plus faces */ + P_s[idim][5] = gownership[2*idim+1]; /* Aufpunkt, plus faces */ + } + + + /* norm. vector surface, yz-minus surface (face0) */ + n_s[0][0]=-1; + n_s[1][0]= 0; + n_s[2][0]= 0; + + /* norm. vector surface, xz-minus surface (face2) */ + n_s[0][2]= 0; + n_s[1][2]=-1; + n_s[2][2]= 0; + + /* norm. vector surface, xy-minus surface (face4) */ + n_s[0][4]= 0; + n_s[1][4]= 0; + n_s[2][4]=-1; + + /* plus surfaces (1,3,5) opposite: (-1)*minus_surface */ + for(idim=0;idim<MAX_DIM;idim++) + { + n_s[idim][1]=(-1)*n_s[idim][0]; + n_s[idim][3]=(-1)*n_s[idim][2]; + n_s[idim][5]=(-1)*n_s[idim][4]; + } + + + printf("P_s (face 0): %d %d %d \n", + P_s[0][0],P_s[1][0],P_s[2][0]); + printf("P_s (face 1): %d %d %d \n\n", + P_s[0][1],P_s[1][1],P_s[2][1]); + + /* Lamda for each surface */ + for (iface=0;iface<2*vardim;iface++) + { + tmp1=0.0; tmp2=0.0; + + for (idim=0;idim<vardim;idim++) + { + tmp1+= (S_l[idim]-P_s[idim][iface])*n_s[idim][iface]; + tmp2+= (u_l[idim]*n_s[idim][iface]); + } + if (tmp2==0) + { + lmd_s[iface]=BAD; + notcoll[iface]=0; + } + else + { + notcoll[iface]=1; + rtmp1 = (CCTK_REAL)(-1.0)*(tmp1)/(tmp2); + if (rtmp1<0.0) lmd_s[iface]=BAD; + else lmd_s[iface]=rtmp1; + printf("lmd_s (face%d) %f \n",iface,rtmp1); + } + } + + /* Find 2 valid intersection with the proc surfaces */ + /* tricky: duplicate points corners and the startpoint */ + /* We have to restrict the line to positive lambda, otherwise + it is getting extremly ugly in finding out the local + points */ + + lmd[0]=BAD; + lmd[1]=BAD; + isec = 0; + + printf("testing for local Startpt: %d %d %d", + S_l[0],S_l[1],S_l[2]); + + if (IsPointOnProc_INT(vardim,gownership,S_l)==1) + { + printf(" ... YES ... valid first\n"); + lmd[0] = 0; + isec = 1; + secpnt[0][0]= S_l[0]; + secpnt[1][0]= S_l[1]; + secpnt[2][0]= S_l[2]; + + } + else + printf(" ... NO \n"); + + + for (iface=0;iface<MAX_FACE; iface++) + { + if (notcoll[iface]==1 && lmd_s[iface]!=BAD) + { + for (idim=0;idim<vardim;idim++) + t_pnt1[idim]=S_l[idim]+lmd_s[iface]*u_l[idim]; + + printf("testing for local intersec. (face %d lamda %+4.2f): %4.2f %4.2f %4.2f\n", + iface, lmd_s[iface], t_pnt1[0],t_pnt1[1],t_pnt1[2]); + + /* check if point is local */ + if (IsPointOnProc_REAL(vardim,gownership,t_pnt1)==1) + { + printf(" ... YES ... "); + + /* first valid face: this point is good */ + if (isec==0) + { + printf(" valid first %d\n",isec); + lmd[0] =lmd_s[iface]; + secpnt[0][0]=t_pnt1[0]; + secpnt[1][0]=t_pnt1[1]; + secpnt[2][0]=t_pnt1[2]; + isec++; + } + /* we already have intersec. check for duplicates (make sure + that we allow the endpoint to be located on the startpoint lmd=0 )*/ + else if ((secpnt[0][0]==t_pnt1[0]) && + (secpnt[1][0]==t_pnt1[1]) && + (secpnt[2][0]==t_pnt1[2])) + { + printf(" but duplicate, skip \n"); + } + else + { + printf(" valid second %d\n",isec); + lmd[1] =lmd_s[iface]; + secpnt[0][1]=t_pnt1[0]; + secpnt[1][1]=t_pnt1[1]; + secpnt[2][1]=t_pnt1[2]; + isec++; + } + } + else + { + printf(" ... NO \n"); + } + } + else printf(" ... skipping face: %d %4.2f \n",iface,lmd_s[iface]); + } + + /* flag if a proc has intersections, order the intersections + if we have any */ + if (lmd[0]!=BAD && lmd[1]!=BAD) + { + prochasline=1; + + /* switch: lmd[0] has to be the closer one to S_l */ + if (lmd[0]>lmd[1]) + { + t_lmd = lmd[0]; + lmd[0] = lmd[1]; + lmd[1] = t_lmd; + for (idim=0;idim<vardim;idim++) + { + t_pnt1[idim] = secpnt[idim][0]; + secpnt[idim][0]= secpnt[idim][1]; + secpnt[idim][1]= t_pnt1[idim]; + } + } + lmd[0]=ceil (lmd[0]); + lmd[1]=floor(lmd[1]); + printf("LOCAL INTERSECTIONS: (ordered in lmd) %4.2f %4.2f \n",lmd[0],lmd[1]); + } + else if (lmd[0]==BAD && lmd[1]==BAD) + { + prochasline = 0; + printf("NO INTERSECTION ON THIS PROC\n"); + } + else if ((lmd[0]==BAD && lmd[1]!=BAD) || + (lmd[1]==BAD && lmd[0]!=BAD)) + { + printf("ERROR: only one local intersection: that is wrong %4.2f %4.2f\n", + lmd[0],lmd[1]); + } + + /* patch if lambdas are greater line length */ + if (len_l>=0) + { + if (len_l>lmd[0] && len_l<lmd[1]) + lmd[1]=floor(len_l); + else if (len_l<lmd[0]) + { + lmd[0]=BAD; + lmd[1]=BAD; + prochasline=0; + } + } + + /* order intersections in lambda-distance from aufpunkt */ + if (prochasline) + { + nlocpoints[0]= ABS(lmd[0]-lmd[1])+1; + for (idim=0;idim<vardim;idim++) + { + locstart[idim] = S_l[idim]+lmd[0]*u_l[idim]; + locend[idim] = S_l[idim]+lmd[1]*u_l[idim]; + } + } + else + { + nlocpoints[0] = 0; + for (idim=0;idim<vardim;idim++) + { + locstart[idim] = BAD; + locend[idim] = BAD; + } + } + + printf("\n>>>>> GLOBAL START/END: (%d %d %d) -> (%d %d %d) = %d\n", + locstart[0],locstart[1],locstart[2], + locend[0],locend[1],locend[2],nlocpoints[0]); + + return(0); + +} diff --git a/src/TestSlab.c b/src/TestSlab.c new file mode 100644 index 0000000..9c6329f --- /dev/null +++ b/src/TestSlab.c @@ -0,0 +1,103 @@ + +#include <stdio.h> +#include <stdlib.h> +#include <math.h> + +#include "cctk.h" +#include "cctk_Parameters.h" +#include "cctk_Arguments.h" +#include "CactusPUGH/PUGH/src/include/pugh.h" +#include "CactusBase/IOUtil/src/ioGH.h" + +#define MAX_DIM 3 +#define MAX_FACE 6 +#define EPS 1e-10 +#define BAD -42 + +#define ABS(a) ((a)<0 ? -(a) : (a)) + +int IsPointOnProc(int dim, int *gownership, int *pnt); +int GetLocalLine(cGH *GH, + int vardim, int *P_l, int *u_l, int *ds, int len_l, + int *locstart, int *locend, int *nlocpoints); + +int CollectLocalData1D(cGH *GH, int vindex, int vtype, int vtypesize, int vtimelvl, + int vdim, int *S_l, int *u_l, int ds, int len_l, + void **loc_dptr, int *dsize); + +void TestSlab(CCTK_ARGUMENTS) { + DECLARE_CCTK_ARGUMENTS + + int gs[3], di[3], ds[1]; /* global start/direction/downsampling*/ + int ls[3], le[3]; /* local start/end */ + int nlocpoints; /* number of points */ + int nprocs,iproc; /* number of procs */ + int *c_nlocpoints; /* collected number of points per proc (size nproc) */ + pGH *pughGH; + int gindex; + int vtype, vindex, vtypesize, vdim, vtimelvl; + void *loc_dptr; + + vindex = CCTK_VarIndex("grid::x"); + gindex = CCTK_GroupIndexFromVarI(vindex); + vtype = CCTK_GroupTypeFromVarI(vindex); + vtypesize= CCTK_VarTypeSize(vtype); + vtimelvl = CCTK_NumTimeLevelsFromVarI (vindex) - 1; + vdim = CCTK_GroupDimI(gindex); + + + gs[0]=1; + gs[1]=5; + gs[2]=9; + + di[0]= 0; + di[1]=-1; + di[2]=-2; + + ds[0]=1; + + CollectLocalData1D(cctkGH, vindex, vtype, vtypesize, vtimelvl, + vdim, gs, di, ds, -1, &loc_dptr, &nlocpoints); + + +} + +void TestGetLocalLine(CCTK_ARGUMENTS) { + DECLARE_CCTK_ARGUMENTS + + int gs[3], di[3], ds[1]; /* global start/direction/downsampling*/ + int ls[3], le[3]; /* local start/end */ + int nlocpoints; /* number of points */ + int nprocs,iproc; /* number of procs */ + int *c_nlocpoints; /* collected number of points per proc (size nproc) */ + pGH *pughGH; + + pughGH = PUGH_pGH (cctkGH); + + gs[0]=1; + gs[1]=5; + gs[2]=9; + + di[0]= 0; + di[1]=-1; + di[2]= 2; + + ds[0]=1; + + GetLocalLine(cctkGH, 3, gs, di, ds, 4, ls, le, &nlocpoints); + printf("TEST: local %d %d %d -> %d %d %d TOTAL: %d\n", + ls[0], ls[1], ls[2], le[0], le[1], le[2], nlocpoints ); + + nprocs = CCTK_nProcs(cctkGH); + c_nlocpoints = (int*) malloc (nprocs*sizeof(int)); + +#ifdef CCTK_MPI + CACTUS_MPI_ERROR(MPI_Gather(&nlocpoints, 1, PUGH_MPI_INT, c_nlocpoints, 1, + PUGH_MPI_INT, 0, pughGH->PUGH_COMM_WORLD)); +#endif + + for (iproc=0;iproc<nprocs;iproc++) + printf(" Proc %d : %d \n",iproc,c_nlocpoints[iproc]); + +} + diff --git a/src/make.code.defn b/src/make.code.defn new file mode 100644 index 0000000..f7f22ae --- /dev/null +++ b/src/make.code.defn @@ -0,0 +1,9 @@ +# Main make.code.defn file for thorn Hyperslab +# $Header$ + +# Source files in this directory +SRCS = GetLocalLine.c TestSlab.c CollectData1D.c + +# Subdirectories containing source files +SUBDIRS = + |