Actual source code: ex28.c
1: static char help[] = "Compare parallel partitioning strategies using matrix graphs\n\n";
3: #include <petscmat.h>
5: int main(int argc, char **args)
6: {
7: MatPartitioning part;
8: IS partis;
9: Mat A = NULL;
10: PetscInt max = -1;
11: PetscInt min = -1;
12: PetscReal balance = 0.0;
13: const PetscInt *ranges = NULL;
14: char filein[PETSC_MAX_PATH_LEN];
15: MPI_Comm comm;
16: PetscMPIInt size;
17: PetscInt p;
18: PetscBool flg;
20: /*load matrix*/
22: PetscInitialize(&argc, &args, (char *)0, help);
23: comm = PETSC_COMM_WORLD;
24: MPI_Comm_size(comm, &size);
25: PetscOptionsGetString(NULL, NULL, "-fin", filein, sizeof(filein), &flg);
26: if (flg) {
27: PetscViewer view;
28: PetscViewerBinaryOpen(comm, filein, FILE_MODE_READ, &view);
29: MatCreate(comm, &A);
30: MatLoad(A, view);
31: PetscViewerDestroy(&view);
32: }
34: /*partition matrix*/
35: MatPartitioningCreate(comm, &part);
36: MatPartitioningSetAdjacency(part, A);
37: MatPartitioningSetFromOptions(part);
38: MatPartitioningApply(part, &partis);
39: MatGetOwnershipRanges(A, &ranges);
40: MatGetSize(A, &min, NULL);
41: for (p = 0; p < size; ++p) {
42: const PetscInt partsize = ranges[p + 1] - ranges[p];
44: max = PetscMax(max, partsize);
45: min = PetscMin(min, partsize);
46: }
47: balance = ((PetscReal)max) / min;
48: PetscPrintf(comm, "ranges: ");
49: for (p = 0; p <= size; ++p) {
50: if (p > 0) PetscPrintf(comm, ", ");
51: PetscPrintf(comm, "%" PetscInt_FMT, ranges[p]);
52: }
53: PetscPrintf(comm, "\n");
54: PetscPrintf(comm, "max:%.0lf min:%.0lf balance:%.11lf\n", (double)max, (double)min, (double)balance);
55: PetscObjectViewFromOptions((PetscObject)partis, NULL, "-partition_view");
56: MatPartitioningDestroy(&part);
57: ISDestroy(&partis);
58: MatDestroy(&A);
59: PetscFinalize();
60: return 0;
61: }