1*c4762a1bSJed Brown static const char help[] = "Test ParMETIS handling of negative weights.\n\n"; 2*c4762a1bSJed Brown 3*c4762a1bSJed Brown /* Test contributed by John Fettig */ 4*c4762a1bSJed Brown 5*c4762a1bSJed Brown /* 6*c4762a1bSJed Brown * This file implements two tests for a bug reported in ParMETIS. These tests are not expected to pass without the 7*c4762a1bSJed Brown * patches in the PETSc distribution of ParMetis. See parmetis.py 8*c4762a1bSJed Brown * 9*c4762a1bSJed Brown * 10*c4762a1bSJed Brown * The bug was reported upstream, but has received no action so far. 11*c4762a1bSJed Brown * 12*c4762a1bSJed Brown * http://glaros.dtc.umn.edu/gkhome/node/837 13*c4762a1bSJed Brown * 14*c4762a1bSJed Brown */ 15*c4762a1bSJed Brown 16*c4762a1bSJed Brown #include <petscsys.h> 17*c4762a1bSJed Brown #include <parmetis.h> 18*c4762a1bSJed Brown 19*c4762a1bSJed Brown #define CHKERRQPARMETIS(n) \ 20*c4762a1bSJed Brown if (n == METIS_ERROR_INPUT) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_LIB,"ParMETIS error due to wrong inputs and/or options"); \ 21*c4762a1bSJed Brown else if (n == METIS_ERROR_MEMORY) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_LIB,"ParMETIS error due to insufficient memory"); \ 22*c4762a1bSJed Brown else if (n == METIS_ERROR) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_LIB,"ParMETIS general error"); \ 23*c4762a1bSJed Brown 24*c4762a1bSJed Brown int main(int argc, char *argv[]) 25*c4762a1bSJed Brown { 26*c4762a1bSJed Brown PetscErrorCode ierr; 27*c4762a1bSJed Brown PetscBool flg; 28*c4762a1bSJed Brown PetscMPIInt rank, size; 29*c4762a1bSJed Brown int i, status; 30*c4762a1bSJed Brown idx_t ni,isize,*vtxdist, *xadj, *adjncy, *vwgt, *part; 31*c4762a1bSJed Brown idx_t wgtflag=0, numflag=0, ncon=1, ndims=3, edgecut=0; 32*c4762a1bSJed Brown idx_t options[5]; 33*c4762a1bSJed Brown PetscReal *xyz; 34*c4762a1bSJed Brown real_t *sxyz, *tpwgts, ubvec[1]; 35*c4762a1bSJed Brown MPI_Comm comm; 36*c4762a1bSJed Brown FILE *fp; 37*c4762a1bSJed Brown char fname[PETSC_MAX_PATH_LEN],prefix[PETSC_MAX_PATH_LEN] = ""; 38*c4762a1bSJed Brown size_t red; 39*c4762a1bSJed Brown 40*c4762a1bSJed Brown ierr = PetscInitialize(&argc,&argv,NULL,help);if (ierr) return ierr; 41*c4762a1bSJed Brown #if defined(PETSC_USE_64BIT_INDICES) 42*c4762a1bSJed Brown ierr = PetscPrintf(PETSC_COMM_WORLD,"This example only works with 32 bit indices\n"); 43*c4762a1bSJed Brown ierr = PetscFinalize(); 44*c4762a1bSJed Brown return ierr; 45*c4762a1bSJed Brown #endif 46*c4762a1bSJed Brown MPI_Comm_rank(PETSC_COMM_WORLD,&rank); 47*c4762a1bSJed Brown MPI_Comm_size(PETSC_COMM_WORLD,&size); 48*c4762a1bSJed Brown 49*c4762a1bSJed Brown ierr = PetscOptionsBegin(PETSC_COMM_WORLD,NULL,"Parmetis test options","");CHKERRQ(ierr); 50*c4762a1bSJed Brown ierr = PetscOptionsString("-prefix","Path and prefix of test file","",prefix,prefix,sizeof(prefix),&flg);CHKERRQ(ierr); 51*c4762a1bSJed Brown if (!flg) SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_USER,"Must specify -prefix"); 52*c4762a1bSJed Brown ierr = PetscOptionsEnd();CHKERRQ(ierr); 53*c4762a1bSJed Brown 54*c4762a1bSJed Brown ierr = PetscMalloc1(size+1,&vtxdist);CHKERRQ(ierr); 55*c4762a1bSJed Brown 56*c4762a1bSJed Brown ierr = PetscSNPrintf(fname,sizeof(fname),"%s.%d.graph",prefix,rank);CHKERRQ(ierr); 57*c4762a1bSJed Brown 58*c4762a1bSJed Brown ierr = PetscFOpen(PETSC_COMM_SELF,fname,"r",&fp);CHKERRQ(ierr); 59*c4762a1bSJed Brown 60*c4762a1bSJed Brown red = fread(vtxdist, sizeof(idx_t), size+1, fp);if (red != (size_t) (size+1)) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SYS,"Unable to read from data file"); 61*c4762a1bSJed Brown 62*c4762a1bSJed Brown ni = vtxdist[rank+1]-vtxdist[rank]; 63*c4762a1bSJed Brown 64*c4762a1bSJed Brown ierr = PetscMalloc1(ni+1,&xadj);CHKERRQ(ierr); 65*c4762a1bSJed Brown 66*c4762a1bSJed Brown red = fread(xadj, sizeof(idx_t), ni+1, fp);if (red != (size_t) (ni+1)) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SYS,"Unable to read from data file"); 67*c4762a1bSJed Brown 68*c4762a1bSJed Brown ierr = PetscMalloc1(xadj[ni],&adjncy);CHKERRQ(ierr); 69*c4762a1bSJed Brown 70*c4762a1bSJed Brown for (i=0; i<ni; i++) { 71*c4762a1bSJed Brown red = fread(&adjncy[xadj[i]], sizeof(idx_t), xadj[i+1]-xadj[i], fp);if (red != (size_t) (xadj[i+1]-xadj[i])) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SYS,"Unable to read from data file"); 72*c4762a1bSJed Brown } 73*c4762a1bSJed Brown 74*c4762a1bSJed Brown ierr = PetscFClose(PETSC_COMM_SELF,fp);CHKERRQ(ierr); 75*c4762a1bSJed Brown 76*c4762a1bSJed Brown ierr = PetscSNPrintf(fname,sizeof(fname),"%s.%d.graph.xyz",prefix,rank);CHKERRQ(ierr); 77*c4762a1bSJed Brown ierr = PetscFOpen(PETSC_COMM_SELF,fname,"r",&fp);CHKERRQ(ierr); 78*c4762a1bSJed Brown 79*c4762a1bSJed Brown ierr = PetscMalloc3(ni*ndims,&xyz,ni,&part,size,&tpwgts);CHKERRQ(ierr); 80*c4762a1bSJed Brown ierr = PetscMalloc1(ni*ndims,&sxyz);CHKERRQ(ierr); 81*c4762a1bSJed Brown 82*c4762a1bSJed Brown red = fread(xyz, sizeof(PetscReal), ndims*ni, fp);if (red != (size_t) (ndims*ni)) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SYS,"Unable to read from data file"); 83*c4762a1bSJed Brown for (i=0; i<ni*ndims; i++) sxyz[i] = (size_t) xyz[i]; 84*c4762a1bSJed Brown 85*c4762a1bSJed Brown ierr = PetscFClose(PETSC_COMM_SELF,fp);CHKERRQ(ierr); 86*c4762a1bSJed Brown 87*c4762a1bSJed Brown vwgt = NULL; 88*c4762a1bSJed Brown 89*c4762a1bSJed Brown for (i = 0; i < size; i++) tpwgts[i] = 1. / size; 90*c4762a1bSJed Brown isize = size; 91*c4762a1bSJed Brown 92*c4762a1bSJed Brown ubvec[0] = 1.05; 93*c4762a1bSJed Brown options[0] = 0; 94*c4762a1bSJed Brown options[1] = 2; 95*c4762a1bSJed Brown options[2] = 15; 96*c4762a1bSJed Brown options[3] = 0; 97*c4762a1bSJed Brown options[4] = 0; 98*c4762a1bSJed Brown 99*c4762a1bSJed Brown ierr = MPI_Comm_dup(MPI_COMM_WORLD, &comm);CHKERRQ(ierr); 100*c4762a1bSJed Brown status = ParMETIS_V3_PartGeomKway(vtxdist, xadj, adjncy, vwgt, NULL, &wgtflag, &numflag, &ndims, sxyz, &ncon, &isize, tpwgts, ubvec,options, &edgecut, part, &comm);CHKERRQPARMETIS(status); 101*c4762a1bSJed Brown ierr = MPI_Comm_free(&comm);CHKERRQ(ierr); 102*c4762a1bSJed Brown 103*c4762a1bSJed Brown ierr = PetscFree(vtxdist);CHKERRQ(ierr); 104*c4762a1bSJed Brown ierr = PetscFree(xadj);CHKERRQ(ierr); 105*c4762a1bSJed Brown ierr = PetscFree(adjncy);CHKERRQ(ierr); 106*c4762a1bSJed Brown ierr = PetscFree3(xyz,part,tpwgts);CHKERRQ(ierr); 107*c4762a1bSJed Brown ierr = PetscFree(sxyz);CHKERRQ(ierr); 108*c4762a1bSJed Brown ierr = PetscFinalize(); 109*c4762a1bSJed Brown return ierr; 110*c4762a1bSJed Brown } 111*c4762a1bSJed Brown 112*c4762a1bSJed Brown 113*c4762a1bSJed Brown /*TEST 114*c4762a1bSJed Brown 115*c4762a1bSJed Brown build: 116*c4762a1bSJed Brown requires: parmetis 117*c4762a1bSJed Brown 118*c4762a1bSJed Brown test: 119*c4762a1bSJed Brown nsize: 2 120*c4762a1bSJed Brown requires: parmetis datafilespath !complex double !define(PETSC_USE_64BIT_INDICES) 121*c4762a1bSJed Brown args: -prefix ${DATAFILESPATH}/parmetis-test/testnp2 122*c4762a1bSJed Brown 123*c4762a1bSJed Brown test: 124*c4762a1bSJed Brown suffix: 2 125*c4762a1bSJed Brown nsize: 4 126*c4762a1bSJed Brown requires: parmetis datafilespath !complex double !define(PETSC_USE_64BIT_INDICES) 127*c4762a1bSJed Brown args: -prefix ${DATAFILESPATH}/parmetis-test/testnp4 128*c4762a1bSJed Brown 129*c4762a1bSJed Brown TEST*/ 130