Actual source code: scotch.c
petsc-3.3-p7 2013-05-11
2: #include <../src/mat/impls/adj/mpi/mpiadj.h> /*I "petscmat.h" I*/
4: EXTERN_C_BEGIN
5: #include <ptscotch.h>
6: EXTERN_C_END
8: typedef struct {
9: double imbalance;
10: SCOTCH_Num strategy;
11: } MatPartitioning_PTScotch;
15: /*@
16: MatPartitioningPTScotchSetImbalance - Sets the value of the load imbalance
17: ratio to be used during strategy selection.
19: Collective on MatPartitioning
21: Input Parameters:
22: + part - the partitioning context
23: - imb - the load imbalance ratio
25: Options Database:
26: . -mat_partitioning_ptscotch_imbalance <imb>
28: Note:
29: Must be in the range [0,1]. The default value is 0.01.
31: Level: advanced
33: .seealso: MatPartitioningPTScotchSetStrategy(), MatPartitioningPTScotchGetImbalance()
34: @*/
35: PetscErrorCode MatPartitioningPTScotchSetImbalance(MatPartitioning part,PetscReal imb)
36: {
42: PetscTryMethod(part,"MatPartitioningPTScotchSetImbalance_C",(MatPartitioning,PetscReal),(part,imb));
43: return(0);
44: }
46: EXTERN_C_BEGIN
49: PetscErrorCode MatPartitioningPTScotchSetImbalance_PTScotch(MatPartitioning part,PetscReal imb)
50: {
51: MatPartitioning_PTScotch *scotch = (MatPartitioning_PTScotch*)part->data;
54: if (imb==PETSC_DEFAULT) scotch->imbalance = 0.01;
55: else {
56: if (imb<0.0 || imb>1.0) SETERRQ(((PetscObject)part)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of imb. Must be in range [0,1]");
57: scotch->imbalance = (double)imb;
58: }
59: return(0);
60: }
61: EXTERN_C_END
65: /*@
66: MatPartitioningPTScotchGetImbalance - Gets the value of the load imbalance
67: ratio used during strategy selection.
69: Not Collective
71: Input Parameter:
72: . part - the partitioning context
74: Output Parameter:
75: . imb - the load imbalance ratio
77: Level: advanced
79: .seealso: MatPartitioningPTScotchSetImbalance()
80: @*/
81: PetscErrorCode MatPartitioningPTScotchGetImbalance(MatPartitioning part,PetscReal *imb)
82: {
88: PetscTryMethod(part,"MatPartitioningPTScotchGetImbalance_C",(MatPartitioning,PetscReal*),(part,imb));
89: return(0);
90: }
92: EXTERN_C_BEGIN
95: PetscErrorCode MatPartitioningPTScotchGetImbalance_PTScotch(MatPartitioning part,PetscReal *imb)
96: {
97: MatPartitioning_PTScotch *scotch = (MatPartitioning_PTScotch*)part->data;
100: *imb = scotch->imbalance;
101: return(0);
102: }
103: EXTERN_C_END
107: /*@
108: MatPartitioningPTScotchSetStrategy - Sets the strategy to be used in PTScotch.
110: Collective on MatPartitioning
112: Input Parameters:
113: + part - the partitioning context
114: - strategy - the strategy, one of
115: .vb
116: MP_PTSCOTCH_QUALITY - Prioritize quality over speed
117: MP_PTSCOTCH_SPEED - Prioritize speed over quality
118: MP_PTSCOTCH_BALANCE - Enforce load balance
119: MP_PTSCOTCH_SAFETY - Avoid methods that may fail
120: MP_PTSCOTCH_SCALABILITY - Favor scalability as much as possible
121: .ve
123: Options Database:
124: . -mat_partitioning_ptscotch_strategy [quality,speed,balance,safety,scalability] - strategy
125:
126: Level: advanced
128: Notes:
129: The default is MP_SCOTCH_QUALITY. See the PTScotch documentation for more information.
131: .seealso: MatPartitioningPTScotchSetImbalance(), MatPartitioningPTScotchGetStrategy()
132: @*/
133: PetscErrorCode MatPartitioningPTScotchSetStrategy(MatPartitioning part,MPPTScotchStrategyType strategy)
134: {
140: PetscTryMethod(part,"MatPartitioningPTScotchSetStrategy_C",(MatPartitioning,MPPTScotchStrategyType),(part,strategy));
141: return(0);
142: }
144: EXTERN_C_BEGIN
147: PetscErrorCode MatPartitioningPTScotchSetStrategy_PTScotch(MatPartitioning part,MPPTScotchStrategyType strategy)
148: {
149: MatPartitioning_PTScotch *scotch = (MatPartitioning_PTScotch*)part->data;
152: switch (strategy) {
153: case MP_PTSCOTCH_QUALITY: scotch->strategy = SCOTCH_STRATQUALITY; break;
154: case MP_PTSCOTCH_SPEED: scotch->strategy = SCOTCH_STRATSPEED; break;
155: case MP_PTSCOTCH_BALANCE: scotch->strategy = SCOTCH_STRATBALANCE; break;
156: case MP_PTSCOTCH_SAFETY: scotch->strategy = SCOTCH_STRATSAFETY; break;
157: case MP_PTSCOTCH_SCALABILITY: scotch->strategy = SCOTCH_STRATSCALABILITY; break;
158: }
159: return(0);
160: }
161: EXTERN_C_END
165: /*@
166: MatPartitioningPTScotchGetStrategy - Gets the strategy used in PTScotch.
168: Not Collective
170: Input Parameter:
171: . part - the partitioning context
173: Output Parameter:
174: . strategy - the strategy
175:
176: Level: advanced
178: .seealso: MatPartitioningPTScotchSetStrategy()
179: @*/
180: PetscErrorCode MatPartitioningPTScotchGetStrategy(MatPartitioning part,MPPTScotchStrategyType *strategy)
181: {
187: PetscTryMethod(part,"MatPartitioningPTScotchGetStrategy_C",(MatPartitioning,MPPTScotchStrategyType*),(part,strategy));
188: return(0);
189: }
191: EXTERN_C_BEGIN
194: PetscErrorCode MatPartitioningPTScotchGetStrategy_PTScotch(MatPartitioning part,MPPTScotchStrategyType *strategy)
195: {
196: MatPartitioning_PTScotch *scotch = (MatPartitioning_PTScotch*)part->data;
199: switch (scotch->strategy) {
200: case SCOTCH_STRATQUALITY: *strategy = MP_PTSCOTCH_QUALITY; break;
201: case SCOTCH_STRATSPEED: *strategy = MP_PTSCOTCH_SPEED; break;
202: case SCOTCH_STRATBALANCE: *strategy = MP_PTSCOTCH_BALANCE; break;
203: case SCOTCH_STRATSAFETY: *strategy = MP_PTSCOTCH_SAFETY; break;
204: case SCOTCH_STRATSCALABILITY: *strategy = MP_PTSCOTCH_SCALABILITY; break;
205: }
206: return(0);
207: }
208: EXTERN_C_END
212: PetscErrorCode MatPartitioningView_PTScotch(MatPartitioning part, PetscViewer viewer)
213: {
214: MatPartitioning_PTScotch *scotch = (MatPartitioning_PTScotch*)part->data;
215: PetscErrorCode ierr;
216: PetscBool isascii;
217: const char *str=0;
218:
220: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
221: if (isascii) {
222: switch (scotch->strategy) {
223: case SCOTCH_STRATQUALITY: str = "Prioritize quality over speed"; break;
224: case SCOTCH_STRATSPEED: str = "Prioritize speed over quality"; break;
225: case SCOTCH_STRATBALANCE: str = "Enforce load balance"; break;
226: case SCOTCH_STRATSAFETY: str = "Avoid methods that may fail"; break;
227: case SCOTCH_STRATSCALABILITY: str = "Favor scalability as much as possible"; break;
228: }
229: PetscViewerASCIIPrintf(viewer," Strategy=%s\n",str);
230: PetscViewerASCIIPrintf(viewer," Load imbalance ratio=%g\n",scotch->imbalance);
231: } else {
232: SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_SUP,"Viewer type %s not supported for PTScotch partitioner",((PetscObject)viewer)->type_name);
233: }
234: return(0);
235: }
239: PetscErrorCode MatPartitioningSetFromOptions_PTScotch(MatPartitioning part)
240: {
241: PetscErrorCode ierr;
242: PetscBool flag;
243: PetscReal r;
244: MatPartitioning_PTScotch *scotch = (MatPartitioning_PTScotch*)part->data;
245: MPPTScotchStrategyType strat;
248: PetscOptionsHead("PTScotch partitioning options");
249: PetscOptionsEnum("-mat_partitioning_ptscotch_strategy","Strategy","MatPartitioningPTScotchSetStrategy",MPPTScotchStrategyTypes,(PetscEnum)MP_PTSCOTCH_QUALITY,(PetscEnum*)&strat,&flag);
250: if (flag) { MatPartitioningPTScotchSetStrategy(part,strat); }
251: PetscOptionsReal("-mat_partitioning_ptscotch_imbalance","Load imbalance ratio","MatPartitioningPTScotchSetImbalance",scotch->imbalance,&r,&flag);
252: if (flag) { MatPartitioningPTScotchSetImbalance(part,r); }
253: PetscOptionsTail();
254: return(0);
255: }
259: PetscErrorCode MatPartitioningApply_PTScotch(MatPartitioning part,IS *partitioning)
260: {
261: MatPartitioning_PTScotch *scotch = (MatPartitioning_PTScotch*)part->data;
262: PetscErrorCode ierr;
263: PetscMPIInt rank;
264: Mat mat = part->adj;
265: Mat_MPIAdj *adj = (Mat_MPIAdj*)mat->data;
266: PetscBool flg;
267: PetscInt i,j,wgtflag=0,bs=1,nold;
268: PetscReal *vwgttab,deltval;
269: SCOTCH_Num *locals,*velotab,*veloloctab,*edloloctab,vertlocnbr,edgelocnbr,nparts=part->n;
270: SCOTCH_Arch archdat;
271: SCOTCH_Dgraph grafdat;
272: SCOTCH_Dmapping mappdat;
273: SCOTCH_Strat stradat;
276: MPI_Comm_rank(((PetscObject)part)->comm,&rank);
277: PetscObjectTypeCompare((PetscObject)mat,MATMPIADJ,&flg);
278: if (!flg) {
279: /* bs indicates if the converted matrix is "reduced" from the original and hence the
280: resulting partition results need to be stretched to match the original matrix */
281: nold = mat->rmap->n;
282: MatConvert(mat,MATMPIADJ,MAT_INITIAL_MATRIX,&mat);
283: bs = nold/mat->rmap->n;
284: adj = (Mat_MPIAdj*)mat->data;
285: }
287: PetscMalloc((mat->rmap->n+1)*sizeof(SCOTCH_Num),&locals);
288: PetscMalloc(nparts*sizeof(PetscReal),&vwgttab);
289: PetscMalloc(nparts*sizeof(SCOTCH_Num),&velotab);
290: for (j=0;j<nparts;j++) {
291: if (part->part_weights) vwgttab[j] = part->part_weights[j]*nparts;
292: else vwgttab[j] = 1.0;
293: }
294: for (i=0;i<nparts;i++) {
295: deltval = PetscAbsReal(vwgttab[i]-floor(vwgttab[i]+0.5));
296: if (deltval>0.01) {
297: for (j=0;j<nparts;j++) vwgttab[j] /= deltval;
298: }
299: }
300: for (i=0;i<nparts;i++)
301: velotab[i] = (SCOTCH_Num) (vwgttab[i] + 0.5);
302: PetscFree(vwgttab);
304: SCOTCH_dgraphInit(&grafdat,((PetscObject)part)->comm);
306: vertlocnbr = mat->rmap->range[rank+1] - mat->rmap->range[rank];
307: edgelocnbr = adj->i[vertlocnbr];
308: veloloctab = (!part->vertex_weights && !(wgtflag & 2)) ? part->vertex_weights: PETSC_NULL;
309: edloloctab = (!adj->values && !(wgtflag & 1)) ? adj->values: PETSC_NULL;
311: SCOTCH_dgraphBuild(&grafdat,0,vertlocnbr,vertlocnbr,adj->i,adj->i+1,veloloctab,
312: PETSC_NULL,edgelocnbr,edgelocnbr,adj->j,PETSC_NULL,edloloctab);
314: #if defined(PETSC_USE_DEBUG)
315: SCOTCH_dgraphCheck(&grafdat);
316: #endif
318: SCOTCH_archInit(&archdat);
319: SCOTCH_stratInit(&stradat);
320: SCOTCH_stratDgraphMapBuild(&stradat,scotch->strategy,nparts,nparts,scotch->imbalance);
322: SCOTCH_archCmpltw(&archdat,nparts,velotab);
323: SCOTCH_dgraphMapInit(&grafdat,&mappdat,&archdat,locals);
324: SCOTCH_dgraphMapCompute(&grafdat,&mappdat,&stradat);
326: SCOTCH_dgraphMapExit (&grafdat,&mappdat);
327: SCOTCH_archExit(&archdat);
328: SCOTCH_stratExit(&stradat);
329: SCOTCH_dgraphExit(&grafdat);
330: PetscFree(velotab);
332: if (bs > 1) {
333: PetscInt *newlocals;
334: PetscMalloc(bs*mat->rmap->n*sizeof(PetscInt),&newlocals);
335: for (i=0;i<mat->rmap->n;i++) {
336: for (j=0;j<bs;j++) {
337: newlocals[bs*i+j] = locals[i];
338: }
339: }
340: PetscFree(locals);
341: ISCreateGeneral(((PetscObject)part)->comm,bs*mat->rmap->n,newlocals,PETSC_OWN_POINTER,partitioning);
342: } else {
343: ISCreateGeneral(((PetscObject)part)->comm,mat->rmap->n,locals,PETSC_OWN_POINTER,partitioning);
344: }
346: if (!flg) {
347: MatDestroy(&mat);
348: }
349: return(0);
350: }
354: PetscErrorCode MatPartitioningDestroy_PTScotch(MatPartitioning part)
355: {
356: MatPartitioning_PTScotch *scotch = (MatPartitioning_PTScotch*)part->data;
357: PetscErrorCode ierr;
360: PetscFree(scotch);
361: /* clear composed functions */
362: PetscObjectComposeFunctionDynamic((PetscObject)part,"MatPartitioningPTScotchSetImbalance_C","",PETSC_NULL);
363: PetscObjectComposeFunctionDynamic((PetscObject)part,"MatPartitioningPTScotchGetImbalance_C","",PETSC_NULL);
364: PetscObjectComposeFunctionDynamic((PetscObject)part,"MatPartitioningPTScotchSetStrategy_C","",PETSC_NULL);
365: PetscObjectComposeFunctionDynamic((PetscObject)part,"MatPartitioningPTScotchGetStrategy_C","",PETSC_NULL);
366: return(0);
367: }
369: /*MC
370: MATPARTITIONINGPTSCOTCH - Creates a partitioning context via the external package SCOTCH.
372: Level: beginner
374: Notes: See http://www.labri.fr/perso/pelegrin/scotch/
376: .keywords: Partitioning, create, context
378: .seealso: MatPartitioningSetType(), MatPartitioningType
379: M*/
381: EXTERN_C_BEGIN
384: PetscErrorCode MatPartitioningCreate_PTScotch(MatPartitioning part)
385: {
386: PetscErrorCode ierr;
387: MatPartitioning_PTScotch *scotch;
390: PetscNewLog(part,MatPartitioning_PTScotch,&scotch);
391: part->data = (void*)scotch;
393: scotch->imbalance = 0.01;
394: scotch->strategy = SCOTCH_STRATQUALITY;
396: part->ops->apply = MatPartitioningApply_PTScotch;
397: part->ops->view = MatPartitioningView_PTScotch;
398: part->ops->setfromoptions = MatPartitioningSetFromOptions_PTScotch;
399: part->ops->destroy = MatPartitioningDestroy_PTScotch;
401: PetscObjectComposeFunctionDynamic((PetscObject)part,"MatPartitioningPTScotchSetImbalance_C","MatPartitioningPTScotchSetImbalance_PTScotch",MatPartitioningPTScotchSetImbalance_PTScotch);
402: PetscObjectComposeFunctionDynamic((PetscObject)part,"MatPartitioningPTScotchGetImbalance_C","MatPartitioningPTScotchGetImbalance_PTScotch",MatPartitioningPTScotchGetImbalance_PTScotch);
403: PetscObjectComposeFunctionDynamic((PetscObject)part,"MatPartitioningPTScotchSetStrategy_C","MatPartitioningPTScotchSetStrategy_PTScotch",MatPartitioningPTScotchSetStrategy_PTScotch);
404: PetscObjectComposeFunctionDynamic((PetscObject)part,"MatPartitioningPTScotchGetStrategy_C","MatPartitioningPTScotchGetStrategy_PTScotch",MatPartitioningPTScotchGetStrategy_PTScotch);
405: return(0);
406: }
407: EXTERN_C_END