Actual source code: scotch.c
1: #define PETSCMAT_DLL
3: #include ../src/mat/impls/adj/mpi/mpiadj.h
5: #ifdef PETSC_HAVE_UNISTD_H
6: #include <unistd.h>
7: #endif
9: #ifdef PETSC_HAVE_STDLIB_H
10: #include <stdlib.h>
11: #endif
13: /*
14: Currently using Scotch-3.4
15: */
17: #include "scotch.h"
20: /*************************************
21: * *
22: * Note: *
23: * *
24: * To make scotch compile I *
25: * modified all old mat->m/M into *
26: * mat->rmap->n/N *
27: * *
28: * Hope I was right *
29: * *
30: *************************************/
31: typedef struct {
32: char arch[PETSC_MAX_PATH_LEN];
33: int multilevel;
34: char strategy[30];
35: int global_method; /* global method */
36: int local_method; /* local method */
37: int nbvtxcoarsed; /* number of vertices for the coarse graph */
38: int map; /* to know if we map on archptr or just partionate the graph */
39: char *mesg_log;
40: char host_list[PETSC_MAX_PATH_LEN];
41: } MatPartitioning_Scotch;
43: #define SIZE_LOG 10000 /* size of buffer for msg_log */
47: static PetscErrorCode MatPartitioningApply_Scotch(MatPartitioning part, IS * partitioning)
48: {
50: int *parttab, *locals = PETSC_NULL, rank, i, size;
51: size_t j;
52: Mat mat = part->adj, matMPI, matSeq;
53: int nb_locals = mat->rmap->n;
54: Mat_MPIAdj *adj = (Mat_MPIAdj *) mat->data;
55: MatPartitioning_Scotch *scotch = (MatPartitioning_Scotch *) part->data;
56: PetscTruth flg;
57: #ifdef PETSC_HAVE_UNISTD_H
58: int fd_stdout, fd_pipe[2], count,err;
59: #endif
63: /* check if the matrix is sequential, use MatGetSubMatrices if necessary */
64: MPI_Comm_size(((PetscObject)mat)->comm, &size);
65: PetscTypeCompare((PetscObject) mat, MATMPIADJ, &flg);
66: if (size > 1) {
67: int M, N;
68: IS isrow, iscol;
69: Mat *A;
71: if (flg) {
72: SETERRQ(0, "Distributed matrix format MPIAdj is not supported for sequential partitioners");
73: }
74: PetscPrintf(((PetscObject)part)->comm, "Converting distributed matrix to sequential: this could be a performance loss\n");
76: MatGetSize(mat, &M, &N);
77: ISCreateStride(PETSC_COMM_SELF, M, 0, 1, &isrow);
78: ISCreateStride(PETSC_COMM_SELF, N, 0, 1, &iscol);
79: MatGetSubMatrices(mat, 1, &isrow, &iscol, MAT_INITIAL_MATRIX, &A);
80: matSeq = *A;
81: PetscFree(A);
82: ISDestroy(isrow);
83: ISDestroy(iscol);
84: } else
85: matSeq = mat;
87: /* convert the the matrix to MPIADJ type if necessary */
88: if (!flg) {
89: MatConvert(matSeq, MATMPIADJ, MAT_INITIAL_MATRIX, &matMPI);
90: } else {
91: matMPI = matSeq;
92: }
94: adj = (Mat_MPIAdj *) matMPI->data; /* finaly adj contains adjacency graph */
96: MPI_Comm_rank(((PetscObject)part)->comm, &rank);
98: {
99: /* definition of Scotch library arguments */
100: SCOTCH_Strat stratptr; /* scotch strategy */
101: SCOTCH_Graph grafptr; /* scotch graph */
102: #if defined(DOES_NOT_COMPILE_DUE_TO_BROKEN_INTERFACE)
103: int vertnbr = mat->rmap->N; /* number of vertices in full graph */
104: int *verttab = adj->i; /* start of edge list for each vertex */
105: int *edgetab = adj->j; /* edge list data */
106: int edgenbr = adj->nz; /* number of edges */
107: int *velotab = NULL; /* not used by petsc interface */
108: int *vlbltab = NULL;
109: int *edlotab = NULL;
110: int flagval = 3; /* (cf doc scotch no weight edge & vertices) */
111: #endif
112: int baseval = 0; /* 0 for C array indexing */
113: char strategy[256];
115: PetscMalloc((mat->rmap->N) * sizeof(int), &parttab);
117: /* redirect output to buffer scotch -> mesg_log */
118: #ifdef PETSC_HAVE_UNISTD_H
119: fd_stdout = dup(1);
120: pipe(fd_pipe);
121: close(1);
122: dup2(fd_pipe[1], 1);
123: PetscMalloc(SIZE_LOG * sizeof(char), &(scotch->mesg_log));
124: #endif
126: /* library call */
128: /* Construction of the scotch graph object */
129: SCOTCH_graphInit(&grafptr);
130: #if defined(DOES_NOT_COMPILE_DUE_TO_BROKEN_INTERFACE)
131: SCOTCH_graphBuild((SCOTCH_Graph *) &grafptr,
132: (const SCOTCH_Num) vertnbr,
133: (const SCOTCH_Num) verttab,
134: (const SCOTCH_Num *)velotab,
135: (const SCOTCH_Num *)vlbltab,
136: (const SCOTCH_Num *)edgenbr,
137: (const SCOTCH_Num *)edgetab,
138: (const SCOTCH_Num) edlotab,
139: (const SCOTCH_Num *)baseval,
140: (const SCOTCH_Num *)flagval);
141: #else
142: SETERRQ(PETSC_ERR_SUP,"Scotch interface currently broken");
143: #endif
144: SCOTCH_graphCheck(&grafptr);
146: /* Construction of the strategy */
147: if (scotch->strategy[0] != 0) {
148: PetscStrcpy(strategy, scotch->strategy);
149: } else {
150: PetscStrcpy(strategy, "b{strat=");
152: if (scotch->multilevel) {
153: /* PetscStrcat(strategy,"m{vert=");
154: sprintf(strategy+strlen(strategy),"%d",scotch->nbvtxcoarsed);
155: PetscStrcat(strategy,",asc="); */
156: sprintf(strategy, "b{strat=m{vert=%d,asc=",
157: scotch->nbvtxcoarsed);
158: } else
159: PetscStrcpy(strategy, "b{strat=");
161: switch (scotch->global_method) {
162: case MP_SCOTCH_GREEDY:
163: PetscStrcat(strategy, "h");
164: break;
165: case MP_SCOTCH_GPS:
166: PetscStrcat(strategy, "g");
167: break;
168: case MP_SCOTCH_GR_GPS:
169: PetscStrcat(strategy, "g|h");
170: }
172: switch (scotch->local_method) {
173: case MP_SCOTCH_KERNIGHAN_LIN:
174: if (scotch->multilevel)
175: PetscStrcat(strategy, ",low=f}");
176: else
177: PetscStrcat(strategy, " f");
178: break;
179: case MP_SCOTCH_NONE:
180: if (scotch->multilevel)
181: PetscStrcat(strategy, ",asc=x}");
182: default:
183: break;
184: }
186: PetscStrcat(strategy, " x}");
187: }
189: PetscPrintf(((PetscObject)part)->comm, "strategy=[%s]\n", strategy);
191: SCOTCH_stratInit(&stratptr);
192: /*
194: TODO: Correct this part
196: Commented because this doesn't exists anymore
198:
199: SCOTCH_stratMap(&stratptr, strategy);
200: */
201: /* check for option mapping */
202: if (!scotch->map) {
203: /* ********************************************
204: * *
205: * TODO: Correct this part *
206: * *
207: * Won't work with this tmp SCOTCH_Strat... *
208: * *
209: * I just modified it to make scotch compile, *
210: * to be able to use PaStiX... *
211: * *
212: **********************************************/
213: #if defined (DOES_NOT_COMPILE_DUE_TO_BROKEN_INTERFACE)
214: SCOTCH_Strat tmp;
215: SCOTCH_graphPart((const SCOTCH_Graph *)&grafptr,
216: (const SCOTCH_Num) &stratptr,
217: (const SCOTCH_Strat *)&tmp, /* The Argument changed from scotch 3.04 it was part->n, */
218: (SCOTCH_Num *) parttab);
219: #else
220: SETERRQ(PETSC_ERR_SUP,"Scotch interface currently broken");
221: #endif
222: PetscPrintf(PETSC_COMM_SELF, "Partition simple without mapping\n");
223: } else {
224: SCOTCH_Graph grafarch;
225: SCOTCH_Num *listtab;
226: SCOTCH_Num listnbr = 0;
227: SCOTCH_Arch archptr; /* file in scotch architecture format */
228: SCOTCH_Strat archstrat;
229: int arch_total_size, *parttab_tmp,err;
230: int cpt;
231: char buf[256];
232: FILE *file1, *file2;
233: char host_buf[256];
235: /* generate the graph that represents the arch */
236: file1 = fopen(scotch->arch, "r");
237: if (!file1) SETERRQ1(PETSC_ERR_FILE_OPEN, "Scotch: unable to open architecture file %s", scotch->arch);
239: SCOTCH_graphInit(&grafarch);
240: SCOTCH_graphLoad(&grafarch, file1, baseval, 3);
242: SCOTCH_graphCheck(&grafarch);
243: SCOTCH_graphSize(&grafarch, &arch_total_size, &cpt);
245: err = fclose(file1);
246: if (err) SETERRQ(PETSC_ERR_SYS,"fclose() failed on file");
248: printf("total size = %d\n", arch_total_size);
250: /* generate the list of nodes currently working */
251: PetscGetHostName(host_buf, 256);
252: PetscStrlen(host_buf, &j);
254: file2 = fopen(scotch->host_list, "r");
255: if (!file2) SETERRQ1(PETSC_ERR_FILE_OPEN, "Scotch: unable to open host list file %s", scotch->host_list);
257: i = -1;
258: flg = PETSC_FALSE;
259: while (!feof(file2) && !flg) {
260: i++;
261: fgets(buf, 256, file2);
262: PetscStrncmp(buf, host_buf, j, &flg);
263: }
264: err = fclose(file2);
265: if (err) SETERRQ(PETSC_ERR_SYS,"fclose() failed on file");
266: if (!flg) SETERRQ1(PETSC_ERR_LIB, "Scotch: unable to find '%s' in host list file", host_buf);
268: listnbr = size;
269: PetscMalloc(sizeof(SCOTCH_Num) * listnbr, &listtab);
271: MPI_Allgather(&i, 1, MPI_INT, listtab, 1, MPI_INT, ((PetscObject)part)->comm);
273: printf("listnbr = %d, listtab = ", listnbr);
274: for (i = 0; i < listnbr; i++)
275: printf("%d ", listtab[i]);
277: printf("\n");
278: err = fflush(stdout);
279: if (err) SETERRQ(PETSC_ERR_SYS,"fflush() failed on file");
281: SCOTCH_stratInit(&archstrat);
282: /**************************************************************
283: * *
284: * TODO: Correct this part *
285: * *
286: * Commented because this doesn't exists anymore *
287: * *
288: * SCOTCH_stratBipart(&archstrat, "fx"); *
289: **************************************************************/
290: SCOTCH_archInit(&archptr);
291: SCOTCH_archBuild(&archptr, &grafarch, listnbr, listtab,
292: &archstrat);
294: PetscMalloc((mat->rmap->N) * sizeof(int), &parttab_tmp);
295: /************************************************************************************
296: * *
297: * TODO: Correct this part *
298: * *
299: * Commented because this doesn't exists anymore *
300: * *
301: * SCOTCH_mapInit(&mappptr, &grafptr, &archptr, parttab_tmp); *
302: * *
303: * SCOTCH_mapCompute(&mappptr, &stratptr); *
304: * *
305: * SCOTCH_mapView(&mappptr, stdout); *
306: ************************************************************************************/
307: /* now we have to set in the real parttab at the good place */
308: /* because the ranks order are different than position in */
309: /* the arch graph */
310: for (i = 0; i < mat->rmap->N; i++) {
311: parttab[i] = parttab_tmp[i];
312: }
314: PetscFree(listtab);
315: SCOTCH_archExit(&archptr);
316: /*************************************************
317: * TODO: Correct this part *
318: * *
319: * Commented because this doesn't exists anymore *
320: * SCOTCH_mapExit(&mappptr); *
321: *************************************************/
322: SCOTCH_stratExit(&archstrat);
323: }
325: /* dump to mesg_log... */
326: #ifdef PETSC_HAVE_UNISTD_H
327: err = fflush(stdout);
328: if (err) SETERRQ(PETSC_ERR_SYS,"fflush() failed on stdout");
330: count = read(fd_pipe[0], scotch->mesg_log, (SIZE_LOG - 1) * sizeof(char));
331: if (count < 0)
332: count = 0;
333: scotch->mesg_log[count] = 0;
334: close(1);
335: dup2(fd_stdout, 1);
336: close(fd_stdout);
337: close(fd_pipe[0]);
338: close(fd_pipe[1]);
339: #endif
341: SCOTCH_graphExit(&grafptr);
342: SCOTCH_stratExit(&stratptr);
343: }
345: if (ierr)
346: SETERRQ(PETSC_ERR_LIB, scotch->mesg_log);
348: /* Creation of the index set */
350: MPI_Comm_rank(((PetscObject)part)->comm, &rank);
351: MPI_Comm_size(((PetscObject)part)->comm, &size);
352: nb_locals = mat->rmap->N / size;
353: locals = parttab + rank * nb_locals;
354: if (rank < mat->rmap->N % size) {
355: nb_locals++;
356: locals += rank;
357: } else
358: locals += mat->rmap->N % size;
359: ISCreateGeneral(((PetscObject)part)->comm, nb_locals, locals, partitioning);
361: /* destroying old objects */
362: PetscFree(parttab);
363: if (matSeq != mat) {
364: MatDestroy(matSeq);
365: }
366: if (matMPI != mat) {
367: MatDestroy(matMPI);
368: }
370: return(0);
371: }
376: PetscErrorCode MatPartitioningView_Scotch(MatPartitioning part, PetscViewer viewer)
377: {
378: MatPartitioning_Scotch *scotch = (MatPartitioning_Scotch *) part->data;
379: PetscErrorCode ierr;
380: PetscMPIInt rank;
381: PetscTruth iascii;
382:
384: MPI_Comm_rank(((PetscObject)part)->comm, &rank);
385: PetscTypeCompare((PetscObject) viewer, PETSC_VIEWER_ASCII, &iascii);
386: if (iascii) {
387: if (!rank && scotch->mesg_log) {
388: PetscViewerASCIIPrintf(viewer, "%s\n", scotch->mesg_log);
389: }
390: } else {
391: SETERRQ1(PETSC_ERR_SUP, "Viewer type %s not supported for this Scotch partitioner",((PetscObject)viewer)->type_name);
392: }
393: return(0);
394: }
398: /*@
399: MatPartitioningScotchSetGlobal - Set method for global partitioning.
401: Input Parameter:
402: . part - the partitioning context
403: . method - MP_SCOTCH_GREED, MP_SCOTCH_GIBBS or MP_SCOTCH_GR_GI (the combination of two)
404: Level: advanced
406: @*/
407: PetscErrorCode MatPartitioningScotchSetGlobal(MatPartitioning part,
408: MPScotchGlobalType global)
409: {
410: MatPartitioning_Scotch *scotch = (MatPartitioning_Scotch *) part->data;
414: switch (global) {
415: case MP_SCOTCH_GREEDY:
416: case MP_SCOTCH_GPS:
417: case MP_SCOTCH_GR_GPS:
418: scotch->global_method = global;
419: break;
420: default:
421: SETERRQ(PETSC_ERR_SUP, "Scotch: Unknown or unsupported option");
422: }
424: return(0);
425: }
429: /*@
430: MatPartitioningScotchSetCoarseLevel - Set the coarse level
431:
432: Input Parameter:
433: . part - the partitioning context
434: . level - the coarse level in range [0.0,1.0]
436: Level: advanced
438: @*/
439: PetscErrorCode MatPartitioningScotchSetCoarseLevel(MatPartitioning part, PetscReal level)
440: {
441: MatPartitioning_Scotch *scotch = (MatPartitioning_Scotch *) part->data;
445: if (level < 0 || level > 1.0) {
446: SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,
447: "Scocth: level of coarsening out of range [0.0-1.0]");
448: } else {
449: /* ********************************************
450: * *
451: * TODO: Correct this part *
452: * *
453: * Won't work with this nbvxtcoarsed *
454: * *
455: * I just modified it to make scotch compile, *
456: * to be able to use PaStiX... *
457: * *
458: **********************************************/
459: scotch->nbvtxcoarsed = 0;
460: /* with scotch 3.0.4 it was : scotch->nbvtxcoarsed = (int)(part->adj->N * level); */
461: }
462: if (scotch->nbvtxcoarsed < 20)
463: scotch->nbvtxcoarsed = 20;
465: return(0);
466: }
470: /*@C
471: MatPartitioningScotchSetStrategy - Set the strategy to be used by Scotch.
472: This is an alternative way of specifying the global method, the local
473: method, the coarse level and the multilevel option.
474:
475: Input Parameter:
476: . part - the partitioning context
477: . level - the strategy in Scotch format. Check Scotch documentation.
479: Level: advanced
481: .seealso: MatPartitioningScotchSetGlobal(), MatPartitioningScotchSetLocal(), MatPartitioningScotchSetCoarseLevel(), MatPartitioningScotchSetMultilevel(),
482: @*/
483: PetscErrorCode MatPartitioningScotchSetStrategy(MatPartitioning part, char *strat)
484: {
485: MatPartitioning_Scotch *scotch = (MatPartitioning_Scotch *) part->data;
489: PetscStrcpy(scotch->strategy, strat);
490: return(0);
491: }
496: /*@
497: MatPartitioningScotchSetLocal - Set method for local partitioning.
499: Input Parameter:
500: . part - the partitioning context
501: . method - MP_SCOTCH_KERNIGHAN_LIN or MP_SCOTCH_NONE
503: Level: advanced
505: @*/
506: PetscErrorCode MatPartitioningScotchSetLocal(MatPartitioning part, MPScotchLocalType local)
507: {
508: MatPartitioning_Scotch *scotch = (MatPartitioning_Scotch *) part->data;
512: switch (local) {
513: case MP_SCOTCH_KERNIGHAN_LIN:
514: case MP_SCOTCH_NONE:
515: scotch->local_method = local;
516: break;
517: default:
518: SETERRQ(PETSC_ERR_ARG_CORRUPT, "Scotch: Unknown or unsupported option");
519: }
521: return(0);
522: }
526: /*@C
527: MatPartitioningScotchSetArch - Specify the file that describes the
528: architecture used for mapping. The format of this file is documented in
529: the Scotch manual.
531: Input Parameter:
532: . part - the partitioning context
533: . file - the name of file
534: Level: advanced
536: Note:
537: If the name is not set, then the default "archgraph.src" is used.
539: .seealso: MatPartitioningScotchSetHostList(),MatPartitioningScotchSetMapping()
540: @*/
541: PetscErrorCode MatPartitioningScotchSetArch(MatPartitioning part, const char *filename)
542: {
543: MatPartitioning_Scotch *scotch = (MatPartitioning_Scotch *) part->data;
547: PetscStrcpy(scotch->arch, filename);
549: return(0);
550: }
554: /*@C
555: MatPartitioningScotchSetHostList - Specify host list file for mapping.
557: Input Parameter:
558: . part - the partitioning context
559: . file - the name of file
561: Level: advanced
563: Notes:
564: The file must consist in a list of hostnames (one per line). These hosts
565: are the ones referred to in the architecture file (see
566: MatPartitioningScotchSetArch()): the first host corresponds to index 0,
567: the second one to index 1, and so on.
568:
569: If the name is not set, then the default "host_list" is used.
570:
571: .seealso: MatPartitioningScotchSetArch(), MatPartitioningScotchSetMapping()
572: @*/
573: PetscErrorCode MatPartitioningScotchSetHostList(MatPartitioning part, const char *filename)
574: {
575: MatPartitioning_Scotch *scotch = (MatPartitioning_Scotch *) part->data;
579: PetscStrcpy(scotch->host_list, filename);
581: return(0);
582: }
586: /*@
587: MatPartitioningScotchSetMultilevel - Activates multilevel partitioning.
589: Input Parameter:
590: . part - the partitioning context
592: Level: advanced
594: @*/
595: PetscErrorCode MatPartitioningScotchSetMultilevel(MatPartitioning part)
596: {
597: MatPartitioning_Scotch *scotch = (MatPartitioning_Scotch *) part->data;
601: scotch->multilevel = 1;
603: return(0);
604: }
609: /*@
610: MatPartitioningScotchSetMapping - Activates architecture mapping for the
611: partitioning algorithm. Architecture mapping tries to enhance the quality
612: of partitioning by using network topology information.
614: Input Parameter:
615: . part - the partitioning context
617: Level: advanced
619: .seealso: MatPartitioningScotchSetArch(),MatPartitioningScotchSetHostList()
620: @*/
621: PetscErrorCode MatPartitioningScotchSetMapping(MatPartitioning part)
622: {
623: MatPartitioning_Scotch *scotch = (MatPartitioning_Scotch *) part->data;
627: scotch->map = 1;
629: return(0);
630: }
634: PetscErrorCode MatPartitioningSetFromOptions_Scotch(MatPartitioning part)
635: {
637: PetscTruth flag;
638: char name[PETSC_MAX_PATH_LEN];
639: int i;
640: PetscReal r;
642: const char *global[] = { "greedy", "gps", "gr_gps" };
643: const char *local[] = { "kernighan-lin", "none" };
646: PetscOptionsHead("Set Scotch partitioning options");
648: PetscOptionsEList("-mat_partitioning_scotch_global",
649: "Global method to use", "MatPartitioningScotchSetGlobal", global, 3,
650: global[0], &i, &flag);
651: if (flag) {
652: MatPartitioningScotchSetGlobal(part, (MPScotchGlobalType)i);
653: }
655: PetscOptionsEList("-mat_partitioning_scotch_local",
656: "Local method to use", "MatPartitioningScotchSetLocal", local, 2,
657: local[0], &i, &flag);
658: if (flag) {
659: MatPartitioningScotchSetLocal(part, (MPScotchLocalType)i);
660: }
662: flag = PETSC_FALSE;
663: PetscOptionsTruth("-mat_partitioning_scotch_mapping", "Use mapping","MatPartitioningScotchSetMapping", flag,&flag,PETSC_NULL);
664: if (flag) {
665: MatPartitioningScotchSetMapping(part);
666: }
668: PetscOptionsString("-mat_partitioning_scotch_arch",
669: "architecture file in scotch format", "MatPartitioningScotchSetArch",
670: "archgraph.src", name, PETSC_MAX_PATH_LEN, &flag);
671: if (flag)
672: MatPartitioningScotchSetArch(part, name);
674: PetscOptionsString("-mat_partitioning_scotch_hosts",
675: "host list filename", "MatPartitioningScotchSetHostList",
676: "host_list", name, PETSC_MAX_PATH_LEN, &flag);
677: if (flag)
678: MatPartitioningScotchSetHostList(part, name);
680: PetscOptionsReal("-mat_partitioning_scotch_coarse_level",
681: "coarse level", "MatPartitioningScotchSetCoarseLevel", 0, &r,
682: &flag);
683: if (flag)
684: MatPartitioningScotchSetCoarseLevel(part, r);
686: flag = PETSC_FALSE;
687: PetscOptionsTruth("-mat_partitioning_scotch_mul", "Use coarse level","MatPartitioningScotchSetMultilevel", flag,&flag,PETSC_NULL);
688: if (flag) {
689: MatPartitioningScotchSetMultilevel(part);
690: }
692: PetscOptionsString("-mat_partitioning_scotch_strategy",
693: "Scotch strategy string",
694: "MatPartitioningScotchSetStrategy", "", name, PETSC_MAX_PATH_LEN,
695: &flag);
696: if (flag)
697: MatPartitioningScotchSetStrategy(part, name);
699: PetscOptionsTail();
700: return(0);
701: }
705: PetscErrorCode MatPartitioningDestroy_Scotch(MatPartitioning part)
706: {
707: MatPartitioning_Scotch *scotch = (MatPartitioning_Scotch *) part->data;
708: PetscErrorCode ierr;
711: PetscFree(scotch->mesg_log);
712: PetscFree(scotch);
713: return(0);
714: }
717: /*MC
718: MAT_PARTITIONING_SCOTCH - Creates a partitioning context via the external package SCOTCH.
720: Collective on MPI_Comm
722: Input Parameter:
723: . part - the partitioning context
725: Options Database Keys:
726: + -mat_partitioning_scotch_global <greedy> (one of) greedy gps gr_gps
727: . -mat_partitioning_scotch_local <kernighan-lin> (one of) kernighan-lin none
728: . -mat_partitioning_scotch_mapping: Use mapping (MatPartitioningScotchSetMapping)
729: . -mat_partitioning_scotch_arch <archgraph.src>: architecture file in scotch format (MatPartitioningScotchSetArch)
730: . -mat_partitioning_scotch_hosts <host_list>: host list filename (MatPartitioningScotchSetHostList)
731: . -mat_partitioning_scotch_coarse_level <0>: coarse level (MatPartitioningScotchSetCoarseLevel)
732: . -mat_partitioning_scotch_mul: Use coarse level (MatPartitioningScotchSetMultilevel)
733: - -mat_partitioning_scotch_strategy <>: Scotch strategy string (MatPartitioningScotchSetStrategy)
735: Level: beginner
737: Notes: See http://www.labri.fr/Perso/~pelegrin/scotch/
739: .keywords: Partitioning, create, context
741: .seealso: MatPartitioningSetType(), MatPartitioningType
743: M*/
748: PetscErrorCode MatPartitioningCreate_Scotch(MatPartitioning part)
749: {
751: MatPartitioning_Scotch *scotch;
754: SETERRQ(PETSC_ERR_SUP,"Sorry, the PETSc interface to scotch has not been updated to the latest Scotch version");
755: PetscNewLog(part,MatPartitioning_Scotch, &scotch);
756: part->data = (void*) scotch;
758: scotch->map = 0;
759: scotch->global_method = MP_SCOTCH_GR_GPS;
760: scotch->local_method = MP_SCOTCH_KERNIGHAN_LIN;
761: PetscStrcpy(scotch->arch, "archgraph.src");
762: scotch->nbvtxcoarsed = 200;
763: PetscStrcpy(scotch->strategy, "");
764: scotch->multilevel = 0;
765: scotch->mesg_log = NULL;
767: PetscStrcpy(scotch->host_list, "host_list");
769: part->ops->apply = MatPartitioningApply_Scotch;
770: part->ops->view = MatPartitioningView_Scotch;
771: part->ops->destroy = MatPartitioningDestroy_Scotch;
772: part->ops->setfromoptions = MatPartitioningSetFromOptions_Scotch;
774: return(0);
775: }