Actual source code: jostle.c

 2:  #include src/mat/impls/adj/mpi/mpiadj.h

  4: #ifdef PETSC_HAVE_UNISTD_H
  5: #include <unistd.h>
  6: #endif

  8: #ifdef PETSC_HAVE_STDLIB_H
  9: #include <stdlib.h>
 10: #endif

 12: #include "petscfix.h"


 16: #include "jostle.h"
 17: /* this function is not declared in 'jostle.h' */


 22: typedef struct {
 23:     int output;
 24:     int coarse_seq;
 25:     int nbvtxcoarsed;           /* number of vertices for the coarse graph */
 26:     char *mesg_log;
 27: } MatPartitioning_Jostle;

 29: #define SIZE_LOG 10000          /* size of buffer for msg_log */

 33: static PetscErrorCode MatPartitioningApply_Jostle(MatPartitioning part, IS * partitioning)
 34: {
 36:     int  size, rank, i;
 37:     Mat mat = part->adj, matMPI;
 38:     Mat_MPIAdj *adj = (Mat_MPIAdj *) mat->data;
 39:     MatPartitioning_Jostle *jostle_struct =
 40:         (MatPartitioning_Jostle *) part->data;
 41:     PetscTruth flg;
 42: #ifdef PETSC_HAVE_UNISTD_H
 43:     int fd_stdout, fd_pipe[2], count;
 44: #endif


 48:     /* check that the number of partitions is equal to the number of processors */
 49:     MPI_Comm_rank(mat->comm, &rank);
 50:     MPI_Comm_size(mat->comm, &size);
 51:     if (part->n != size) {
 52:         SETERRQ(PETSC_ERR_SUP, "Supports exactly one domain per processor");
 53:     }

 55:     /* convert adjacency matrix to MPIAdj if needed*/
 56:     PetscTypeCompare((PetscObject) mat, MATMPIADJ, &flg);
 57:     if (!flg) {
 58:         MatConvert(mat, MATMPIADJ, &matMPI);
 59:     } else
 60:         matMPI = mat;

 62:     adj = (Mat_MPIAdj *) matMPI->data;  /* adj contains adjacency graph */
 63:     {
 64:         /* definition of Jostle library arguments */
 65:         int nnodes = matMPI->M; /* number of vertices in full graph */
 66:         int offset = 0;         /* 0 for C array indexing */
 67:         int core = matMPI->m;
 68:         int halo = 0;           /* obsolete with contiguous format */
 69:         int *index_jostle;      /* contribution of each processor */
 70:         int nparts = part->n;
 71:         int *part_wt = NULL;

 73:         int *partition;         /* set number of each vtx (length n) */
 74:         int *degree;            /* degree for each core nodes */
 75:         int *edges = adj->j;
 76:         int *node_wt = NULL;    /* nodes weights */
 77:         int *edge_wt = NULL;    /* edges weights */
 78:         double *coords = NULL;  /* not used (cf jostle documentation) */

 80:         int local_nedges = adj->nz;
 81:         int dimension = 0;      /* not used */
 82:         int output_level = jostle_struct->output;
 83:         char env_str[256];

 85:         /* allocate index_jostle */
 86:         PetscMalloc(nparts * sizeof(int), &index_jostle);

 88:         /* compute number of core nodes for each one */
 89:         for (i = 0; i < nparts - 1; i++)
 90:             index_jostle[i] = adj->rowners[i + 1] - adj->rowners[i];
 91:         index_jostle[nparts - 1] = nnodes - adj->rowners[nparts - 1];

 93:         /* allocate the partition vector */
 94:         PetscMalloc(core * sizeof(int), &partition);

 96:         /* build the degree vector and the local_nedges value */
 97:         PetscMalloc(core * sizeof(int), &degree);
 98:         for (i = 0; i < core; i++)
 99:             degree[i] = adj->i[i + 1] - adj->i[i];

101:         /* library call */
102:         pjostle_init(&size, &rank);
103:         pjostle_comm(&matMPI->comm);
104:         jostle_env("format = contiguous");
105:         jostle_env("timer = off");

107:         sprintf(env_str, "threshold = %d", jostle_struct->nbvtxcoarsed);
108:         jostle_env(env_str);

110:         if (jostle_struct->coarse_seq)
111:             jostle_env("matching = local");

113:         /* redirect output */
114: #ifdef PETSC_HAVE_UNISTD_H
115:         fd_stdout = dup(1);
116:         pipe(fd_pipe);
117:         close(1);
118:         dup2(fd_pipe[1], 1);
119: #endif

121:         pjostle(&nnodes, &offset, &core, &halo, index_jostle, degree, node_wt,
122:             partition, &local_nedges, edges, edge_wt, &nparts,
123:             part_wt, &output_level, &dimension, coords);

125:         printf("Jostle Partitioner statistics\ncut : %d, balance : %f, runtime : %f, mem used : %d\n",
126:             jostle_cut(), jostle_bal(), jostle_tim(), jostle_mem());

128: #ifdef PETSC_HAVE_UNISTD_H
129:         PetscMalloc(SIZE_LOG * sizeof(char), &(jostle_struct->mesg_log));
130:         fflush(stdout);
131:         count = read(fd_pipe[0], jostle_struct->mesg_log, (SIZE_LOG - 1) * sizeof(char));
132:         if (count < 0)
133:             count = 0;
134:         jostle_struct->mesg_log[count] = 0;
135:         close(1);
136:         dup2(fd_stdout, 1);
137:         close(fd_stdout);
138:         close(fd_pipe[0]);
139:         close(fd_pipe[1]);
140: #endif

142:         /* We free the memory used by jostle */
143:         PetscFree(index_jostle);
144:         PetscFree(degree);

146:         /* Creation of the index set */
147:         ISCreateGeneral(part->comm, mat->m, partition, partitioning);

149:         if (matMPI != mat) {
150:             MatDestroy(matMPI);
151:         }

153:         PetscFree(partition);
154:     }

156:     return(0);
157: }


162: PetscErrorCode MatPartitioningView_Jostle(MatPartitioning part, PetscViewer viewer)
163: {
164:   MatPartitioning_Jostle *jostle_struct = (MatPartitioning_Jostle *) part->data;
165:   PetscErrorCode         ierr;
166:   PetscTruth             iascii;

169:   PetscTypeCompare((PetscObject) viewer, PETSC_VIEWER_ASCII, &iascii);
170:   if (iascii) {
171:     if (jostle_struct->mesg_log) {
172:       PetscViewerASCIIPrintf(viewer, "%s\n", jostle_struct->mesg_log);
173:     }
174:   } else {
175:     SETERRQ1(PETSC_ERR_SUP, "Viewer type %s not supported for this Jostle partitioner",((PetscObject) viewer)->type_name);
176:   }
177:   return(0);
178: }

182: /*@
183:     MatPartitioningJostleSetCoarseLevel - Set the coarse level 
184:     
185:   Input Parameter:
186: .  part - the partitioning context
187: .  level - the coarse level in range [0.0,1.0]

189:    Level: advanced

191: @*/
192: PetscErrorCode MatPartitioningJostleSetCoarseLevel(MatPartitioning part, PetscReal level)
193: {
194:     MatPartitioning_Jostle *jostle_struct = (MatPartitioning_Jostle *) part->data;


198:     if (level < 0.0 || level > 1.0) {
199:         SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,
200:             "Jostle: level of coarsening out of range [0.0-1.0]");
201:     } else
202:         jostle_struct->nbvtxcoarsed = (int)(part->adj->N * level);

204:     if (jostle_struct->nbvtxcoarsed < 20)
205:         jostle_struct->nbvtxcoarsed = 20;

207:     return(0);
208: }

212: /*@
213:      MatPartitioningJostleSetCoarseSequential - Use the sequential code to 
214:          do the partitioning of the coarse grid.

216:   Input Parameter:
217: .  part - the partitioning context

219:    Level: advanced

221: @*/
222: PetscErrorCode MatPartitioningJostleSetCoarseSequential(MatPartitioning part)
223: {
224:     MatPartitioning_Jostle *jostle_struct =
225:         (MatPartitioning_Jostle *) part->data;
227:     jostle_struct->coarse_seq = 1;
228:     return(0);
229: }

233: PetscErrorCode MatPartitioningSetFromOptions_Jostle(MatPartitioning part)
234: {
236:     PetscTruth flag;
237:     PetscReal level;

240:     PetscOptionsHead("Set Jostle partitioning options");

242:     PetscOptionsReal("-mat_partitioning_jostle_coarse_level",
243:         "Coarse level", "MatPartitioningJostleSetCoarseLevel", 0,
244:         &level, &flag);
245:     if (flag)
246:         MatPartitioningJostleSetCoarseLevel(part, level);

248:     PetscOptionsName("-mat_partitioning_jostle_coarse_sequential",
249:         "Use sequential coarse partitioner",
250:         "MatPartitioningJostleSetCoarseSequential", &flag);
251:     if (flag) {
252:         MatPartitioningJostleSetCoarseSequential(part);
253:     }

255:     PetscOptionsTail();
256:     return(0);
257: }


262: PetscErrorCode MatPartitioningDestroy_Jostle(MatPartitioning part)
263: {
264:     MatPartitioning_Jostle *jostle_struct =
265:         (MatPartitioning_Jostle *) part->data;


270:     if (jostle_struct->mesg_log) {
271:         PetscFree(jostle_struct->mesg_log);
272:     }
273:     PetscFree(jostle_struct);

275:     return(0);
276: }

281: PetscErrorCode MatPartitioningCreate_Jostle(MatPartitioning part)
282: {
284:     MatPartitioning_Jostle *jostle_struct;

287:     PetscNew(MatPartitioning_Jostle, &jostle_struct);

289:     jostle_struct->nbvtxcoarsed = 20;
290:     jostle_struct->output = 0;
291:     jostle_struct->coarse_seq = 0;
292:     jostle_struct->mesg_log = NULL;

294:     part->ops->apply = MatPartitioningApply_Jostle;
295:     part->ops->view = MatPartitioningView_Jostle;
296:     part->ops->destroy = MatPartitioningDestroy_Jostle;
297:     part->ops->setfromoptions = MatPartitioningSetFromOptions_Jostle;
298:     part->data = (void*) jostle_struct;

300:     return(0);
301: }