Actual source code: scotch.c
petsc-3.8.4 2018-03-24
2: #include <../src/mat/impls/adj/mpi/mpiadj.h>
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;
13: /*@
14: MatPartitioningPTScotchSetImbalance - Sets the value of the load imbalance
15: ratio to be used during strategy selection.
17: Collective on MatPartitioning
19: Input Parameters:
20: + part - the partitioning context
21: - imb - the load imbalance ratio
23: Options Database:
24: . -mat_partitioning_ptscotch_imbalance <imb>
26: Note:
27: Must be in the range [0,1]. The default value is 0.01.
29: Level: advanced
31: .seealso: MatPartitioningPTScotchSetStrategy(), MatPartitioningPTScotchGetImbalance()
32: @*/
33: PetscErrorCode MatPartitioningPTScotchSetImbalance(MatPartitioning part,PetscReal imb)
34: {
40: PetscTryMethod(part,"MatPartitioningPTScotchSetImbalance_C",(MatPartitioning,PetscReal),(part,imb));
41: return(0);
42: }
44: PetscErrorCode MatPartitioningPTScotchSetImbalance_PTScotch(MatPartitioning part,PetscReal imb)
45: {
46: MatPartitioning_PTScotch *scotch = (MatPartitioning_PTScotch*)part->data;
49: if (imb==PETSC_DEFAULT) scotch->imbalance = 0.01;
50: else {
51: if (imb<0.0 || imb>1.0) SETERRQ(PetscObjectComm((PetscObject)part),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of imb. Must be in range [0,1]");
52: scotch->imbalance = (double)imb;
53: }
54: return(0);
55: }
57: /*@
58: MatPartitioningPTScotchGetImbalance - Gets the value of the load imbalance
59: ratio used during strategy selection.
61: Not Collective
63: Input Parameter:
64: . part - the partitioning context
66: Output Parameter:
67: . imb - the load imbalance ratio
69: Level: advanced
71: .seealso: MatPartitioningPTScotchSetImbalance()
72: @*/
73: PetscErrorCode MatPartitioningPTScotchGetImbalance(MatPartitioning part,PetscReal *imb)
74: {
80: PetscUseMethod(part,"MatPartitioningPTScotchGetImbalance_C",(MatPartitioning,PetscReal*),(part,imb));
81: return(0);
82: }
84: PetscErrorCode MatPartitioningPTScotchGetImbalance_PTScotch(MatPartitioning part,PetscReal *imb)
85: {
86: MatPartitioning_PTScotch *scotch = (MatPartitioning_PTScotch*)part->data;
89: *imb = scotch->imbalance;
90: return(0);
91: }
93: /*@
94: MatPartitioningPTScotchSetStrategy - Sets the strategy to be used in PTScotch.
96: Collective on MatPartitioning
98: Input Parameters:
99: + part - the partitioning context
100: - strategy - the strategy, one of
101: .vb
102: MP_PTSCOTCH_QUALITY - Prioritize quality over speed
103: MP_PTSCOTCH_SPEED - Prioritize speed over quality
104: MP_PTSCOTCH_BALANCE - Enforce load balance
105: MP_PTSCOTCH_SAFETY - Avoid methods that may fail
106: MP_PTSCOTCH_SCALABILITY - Favor scalability as much as possible
107: .ve
109: Options Database:
110: . -mat_partitioning_ptscotch_strategy [quality,speed,balance,safety,scalability] - strategy
112: Level: advanced
114: Notes:
115: The default is MP_SCOTCH_QUALITY. See the PTScotch documentation for more information.
117: .seealso: MatPartitioningPTScotchSetImbalance(), MatPartitioningPTScotchGetStrategy()
118: @*/
119: PetscErrorCode MatPartitioningPTScotchSetStrategy(MatPartitioning part,MPPTScotchStrategyType strategy)
120: {
126: PetscTryMethod(part,"MatPartitioningPTScotchSetStrategy_C",(MatPartitioning,MPPTScotchStrategyType),(part,strategy));
127: return(0);
128: }
130: PetscErrorCode MatPartitioningPTScotchSetStrategy_PTScotch(MatPartitioning part,MPPTScotchStrategyType strategy)
131: {
132: MatPartitioning_PTScotch *scotch = (MatPartitioning_PTScotch*)part->data;
135: switch (strategy) {
136: case MP_PTSCOTCH_QUALITY: scotch->strategy = SCOTCH_STRATQUALITY; break;
137: case MP_PTSCOTCH_SPEED: scotch->strategy = SCOTCH_STRATSPEED; break;
138: case MP_PTSCOTCH_BALANCE: scotch->strategy = SCOTCH_STRATBALANCE; break;
139: case MP_PTSCOTCH_SAFETY: scotch->strategy = SCOTCH_STRATSAFETY; break;
140: case MP_PTSCOTCH_SCALABILITY: scotch->strategy = SCOTCH_STRATSCALABILITY; break;
141: }
142: return(0);
143: }
145: /*@
146: MatPartitioningPTScotchGetStrategy - Gets the strategy used in PTScotch.
148: Not Collective
150: Input Parameter:
151: . part - the partitioning context
153: Output Parameter:
154: . strategy - the strategy
156: Level: advanced
158: .seealso: MatPartitioningPTScotchSetStrategy()
159: @*/
160: PetscErrorCode MatPartitioningPTScotchGetStrategy(MatPartitioning part,MPPTScotchStrategyType *strategy)
161: {
167: PetscUseMethod(part,"MatPartitioningPTScotchGetStrategy_C",(MatPartitioning,MPPTScotchStrategyType*),(part,strategy));
168: return(0);
169: }
171: PetscErrorCode MatPartitioningPTScotchGetStrategy_PTScotch(MatPartitioning part,MPPTScotchStrategyType *strategy)
172: {
173: MatPartitioning_PTScotch *scotch = (MatPartitioning_PTScotch*)part->data;
176: switch (scotch->strategy) {
177: case SCOTCH_STRATQUALITY: *strategy = MP_PTSCOTCH_QUALITY; break;
178: case SCOTCH_STRATSPEED: *strategy = MP_PTSCOTCH_SPEED; break;
179: case SCOTCH_STRATBALANCE: *strategy = MP_PTSCOTCH_BALANCE; break;
180: case SCOTCH_STRATSAFETY: *strategy = MP_PTSCOTCH_SAFETY; break;
181: case SCOTCH_STRATSCALABILITY: *strategy = MP_PTSCOTCH_SCALABILITY; break;
182: }
183: return(0);
184: }
186: PetscErrorCode MatPartitioningView_PTScotch(MatPartitioning part, PetscViewer viewer)
187: {
188: MatPartitioning_PTScotch *scotch = (MatPartitioning_PTScotch*)part->data;
189: PetscErrorCode ierr;
190: PetscBool isascii;
191: const char *str=0;
194: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
195: if (isascii) {
196: switch (scotch->strategy) {
197: case SCOTCH_STRATQUALITY: str = "Prioritize quality over speed"; break;
198: case SCOTCH_STRATSPEED: str = "Prioritize speed over quality"; break;
199: case SCOTCH_STRATBALANCE: str = "Enforce load balance"; break;
200: case SCOTCH_STRATSAFETY: str = "Avoid methods that may fail"; break;
201: case SCOTCH_STRATSCALABILITY: str = "Favor scalability as much as possible"; break;
202: }
203: PetscViewerASCIIPrintf(viewer," Strategy=%s\n",str);
204: PetscViewerASCIIPrintf(viewer," Load imbalance ratio=%g\n",scotch->imbalance);
205: }
206: return(0);
207: }
209: PetscErrorCode MatPartitioningSetFromOptions_PTScotch(PetscOptionItems *PetscOptionsObject,MatPartitioning part)
210: {
211: PetscErrorCode ierr;
212: PetscBool flag;
213: PetscReal r;
214: MatPartitioning_PTScotch *scotch = (MatPartitioning_PTScotch*)part->data;
215: MPPTScotchStrategyType strat;
218: MatPartitioningPTScotchGetStrategy(part,&strat);
219: PetscOptionsHead(PetscOptionsObject,"PTScotch partitioning options");
220: PetscOptionsEnum("-mat_partitioning_ptscotch_strategy","Strategy","MatPartitioningPTScotchSetStrategy",MPPTScotchStrategyTypes,(PetscEnum)strat,(PetscEnum*)&strat,&flag);
221: if (flag) { MatPartitioningPTScotchSetStrategy(part,strat); }
222: PetscOptionsReal("-mat_partitioning_ptscotch_imbalance","Load imbalance ratio","MatPartitioningPTScotchSetImbalance",scotch->imbalance,&r,&flag);
223: if (flag) { MatPartitioningPTScotchSetImbalance(part,r); }
224: PetscOptionsTail();
225: return(0);
226: }
228: PetscErrorCode MatPartitioningApply_PTScotch(MatPartitioning part,IS *partitioning)
229: {
230: MatPartitioning_PTScotch *scotch = (MatPartitioning_PTScotch*)part->data;
231: PetscErrorCode ierr;
232: PetscMPIInt rank;
233: Mat mat = part->adj;
234: Mat_MPIAdj *adj = (Mat_MPIAdj*)mat->data;
235: PetscBool flg;
236: PetscInt i,j,wgtflag=0,bs=1,nold;
237: PetscReal *vwgttab,deltval;
238: SCOTCH_Num *locals,*velotab,*veloloctab,*edloloctab,vertlocnbr,edgelocnbr,nparts=part->n;
239: SCOTCH_Arch archdat;
240: SCOTCH_Dgraph grafdat;
241: SCOTCH_Dmapping mappdat;
242: SCOTCH_Strat stradat;
245: MPI_Comm_rank(PetscObjectComm((PetscObject)part),&rank);
246: PetscObjectTypeCompare((PetscObject)mat,MATMPIADJ,&flg);
247: if (!flg) {
248: /* bs indicates if the converted matrix is "reduced" from the original and hence the
249: resulting partition results need to be stretched to match the original matrix */
250: nold = mat->rmap->n;
251: MatConvert(mat,MATMPIADJ,MAT_INITIAL_MATRIX,&mat);
252: if (mat->rmap->n > 0) bs = nold/mat->rmap->n;
253: adj = (Mat_MPIAdj*)mat->data;
254: }
256: PetscMalloc1(mat->rmap->n+1,&locals);
257: PetscMalloc1(nparts,&vwgttab);
258: PetscMalloc1(nparts,&velotab);
259: for (j=0; j<nparts; j++) {
260: if (part->part_weights) vwgttab[j] = part->part_weights[j]*nparts;
261: else vwgttab[j] = 1.0;
262: }
263: for (i=0; i<nparts; i++) {
264: deltval = PetscAbsReal(vwgttab[i]-PetscFloorReal(vwgttab[i]+0.5));
265: if (deltval>0.01) {
266: for (j=0; j<nparts; j++) vwgttab[j] /= deltval;
267: }
268: }
269: for (i=0; i<nparts; i++) velotab[i] = (SCOTCH_Num)(vwgttab[i] + 0.5);
270: PetscFree(vwgttab);
272: SCOTCH_dgraphInit(&grafdat,PetscObjectComm((PetscObject)part));
274: vertlocnbr = mat->rmap->range[rank+1] - mat->rmap->range[rank];
275: edgelocnbr = adj->i[vertlocnbr];
276: veloloctab = (!part->vertex_weights && !(wgtflag & 2)) ? part->vertex_weights : NULL;
277: edloloctab = (!adj->values && !(wgtflag & 1)) ? adj->values : NULL;
279: SCOTCH_dgraphBuild(&grafdat,0,vertlocnbr,vertlocnbr,adj->i,adj->i+1,veloloctab,
280: NULL,edgelocnbr,edgelocnbr,adj->j,NULL,edloloctab);
282: #if defined(PETSC_USE_DEBUG)
283: SCOTCH_dgraphCheck(&grafdat);
284: #endif
286: SCOTCH_archInit(&archdat);
287: SCOTCH_stratInit(&stradat);
288: SCOTCH_stratDgraphMapBuild(&stradat,scotch->strategy,nparts,nparts,scotch->imbalance);
290: SCOTCH_archCmpltw(&archdat,nparts,velotab);
291: SCOTCH_dgraphMapInit(&grafdat,&mappdat,&archdat,locals);
292: SCOTCH_dgraphMapCompute(&grafdat,&mappdat,&stradat);
294: SCOTCH_dgraphMapExit (&grafdat,&mappdat);
295: SCOTCH_archExit(&archdat);
296: SCOTCH_stratExit(&stradat);
297: SCOTCH_dgraphExit(&grafdat);
298: PetscFree(velotab);
300: if (bs > 1) {
301: PetscInt *newlocals;
302: PetscMalloc1(bs*mat->rmap->n,&newlocals);
303: for (i=0;i<mat->rmap->n;i++) {
304: for (j=0;j<bs;j++) {
305: newlocals[bs*i+j] = locals[i];
306: }
307: }
308: PetscFree(locals);
309: ISCreateGeneral(PetscObjectComm((PetscObject)part),bs*mat->rmap->n,newlocals,PETSC_OWN_POINTER,partitioning);
310: } else {
311: ISCreateGeneral(PetscObjectComm((PetscObject)part),mat->rmap->n,locals,PETSC_OWN_POINTER,partitioning);
312: }
314: if (!flg) {
315: MatDestroy(&mat);
316: }
317: return(0);
318: }
320: PetscErrorCode MatPartitioningDestroy_PTScotch(MatPartitioning part)
321: {
322: MatPartitioning_PTScotch *scotch = (MatPartitioning_PTScotch*)part->data;
323: PetscErrorCode ierr;
326: PetscFree(scotch);
327: /* clear composed functions */
328: PetscObjectComposeFunction((PetscObject)part,"MatPartitioningPTScotchSetImbalance_C",NULL);
329: PetscObjectComposeFunction((PetscObject)part,"MatPartitioningPTScotchGetImbalance_C",NULL);
330: PetscObjectComposeFunction((PetscObject)part,"MatPartitioningPTScotchSetStrategy_C",NULL);
331: PetscObjectComposeFunction((PetscObject)part,"MatPartitioningPTScotchGetStrategy_C",NULL);
332: return(0);
333: }
335: /*MC
336: MATPARTITIONINGPTSCOTCH - Creates a partitioning context via the external package SCOTCH.
338: Level: beginner
340: Notes: See http://www.labri.fr/perso/pelegrin/scotch/
342: .keywords: Partitioning, create, context
344: .seealso: MatPartitioningSetType(), MatPartitioningType
345: M*/
347: PETSC_EXTERN PetscErrorCode MatPartitioningCreate_PTScotch(MatPartitioning part)
348: {
349: PetscErrorCode ierr;
350: MatPartitioning_PTScotch *scotch;
353: PetscNewLog(part,&scotch);
354: part->data = (void*)scotch;
356: scotch->imbalance = 0.01;
357: scotch->strategy = SCOTCH_STRATQUALITY;
359: part->ops->apply = MatPartitioningApply_PTScotch;
360: part->ops->view = MatPartitioningView_PTScotch;
361: part->ops->setfromoptions = MatPartitioningSetFromOptions_PTScotch;
362: part->ops->destroy = MatPartitioningDestroy_PTScotch;
364: PetscObjectComposeFunction((PetscObject)part,"MatPartitioningPTScotchSetImbalance_C",MatPartitioningPTScotchSetImbalance_PTScotch);
365: PetscObjectComposeFunction((PetscObject)part,"MatPartitioningPTScotchGetImbalance_C",MatPartitioningPTScotchGetImbalance_PTScotch);
366: PetscObjectComposeFunction((PetscObject)part,"MatPartitioningPTScotchSetStrategy_C",MatPartitioningPTScotchSetStrategy_PTScotch);
367: PetscObjectComposeFunction((PetscObject)part,"MatPartitioningPTScotchGetStrategy_C",MatPartitioningPTScotchGetStrategy_PTScotch);
368: return(0);
369: }