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), °ree);
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: }