1 #include <petsc/private/matimpl.h> /*I "petscmat.h" I*/
2 #include <petscsf.h>
3
MatColoringCreateBipartiteGraph(MatColoring mc,PetscSF * etoc,PetscSF * etor)4 PETSC_EXTERN PetscErrorCode MatColoringCreateBipartiteGraph(MatColoring mc, PetscSF *etoc, PetscSF *etor)
5 {
6 PetscInt nentries, ncolentries, idx;
7 PetscInt i, j, rs, re, cs, ce, cn;
8 PetscInt *rowleaf, *colleaf, *rowdata;
9 PetscInt ncol;
10 const PetscScalar *vcol;
11 const PetscInt *icol;
12 const PetscInt *coldegrees, *rowdegrees;
13 Mat m = mc->mat;
14
15 PetscFunctionBegin;
16 PetscCall(MatGetOwnershipRange(m, &rs, &re));
17 PetscCall(MatGetOwnershipRangeColumn(m, &cs, &ce));
18 cn = ce - cs;
19 nentries = 0;
20 for (i = rs; i < re; i++) {
21 PetscCall(MatGetRow(m, i, &ncol, NULL, &vcol));
22 for (j = 0; j < ncol; j++) nentries++;
23 PetscCall(MatRestoreRow(m, i, &ncol, NULL, &vcol));
24 }
25 PetscCall(PetscMalloc1(nentries, &rowleaf));
26 PetscCall(PetscMalloc1(nentries, &rowdata));
27 idx = 0;
28 for (i = rs; i < re; i++) {
29 PetscCall(MatGetRow(m, i, &ncol, &icol, &vcol));
30 for (j = 0; j < ncol; j++) {
31 rowleaf[idx] = icol[j];
32 rowdata[idx] = i;
33 idx++;
34 }
35 PetscCall(MatRestoreRow(m, i, &ncol, &icol, &vcol));
36 }
37 PetscCheck(idx == nentries, PetscObjectComm((PetscObject)m), PETSC_ERR_NOT_CONVERGED, "Bad number of entries %" PetscInt_FMT " vs %" PetscInt_FMT, idx, nentries);
38 PetscCall(PetscSFCreate(PetscObjectComm((PetscObject)m), etoc));
39 PetscCall(PetscSFCreate(PetscObjectComm((PetscObject)m), etor));
40
41 PetscCall(PetscSFSetGraphLayout(*etoc, m->cmap, nentries, NULL, PETSC_COPY_VALUES, rowleaf));
42 PetscCall(PetscSFSetFromOptions(*etoc));
43
44 /* determine the number of entries in the column matrix */
45 PetscCall(PetscLogEventBegin(MATCOLORING_Comm, *etoc, 0, 0, 0));
46 PetscCall(PetscSFComputeDegreeBegin(*etoc, &coldegrees));
47 PetscCall(PetscSFComputeDegreeEnd(*etoc, &coldegrees));
48 PetscCall(PetscLogEventEnd(MATCOLORING_Comm, *etoc, 0, 0, 0));
49 ncolentries = 0;
50 for (i = 0; i < cn; i++) ncolentries += coldegrees[i];
51 PetscCall(PetscMalloc1(ncolentries, &colleaf));
52
53 /* create the one going the other way by building the leaf set */
54 PetscCall(PetscLogEventBegin(MATCOLORING_Comm, *etoc, 0, 0, 0));
55 PetscCall(PetscSFGatherBegin(*etoc, MPIU_INT, rowdata, colleaf));
56 PetscCall(PetscSFGatherEnd(*etoc, MPIU_INT, rowdata, colleaf));
57 PetscCall(PetscLogEventEnd(MATCOLORING_Comm, *etoc, 0, 0, 0));
58
59 /* this one takes mat entries in *columns* to rows -- you never have to actually be able to order the leaf entries. */
60 PetscCall(PetscSFSetGraphLayout(*etor, m->rmap, ncolentries, NULL, PETSC_COPY_VALUES, colleaf));
61 PetscCall(PetscSFSetFromOptions(*etor));
62
63 PetscCall(PetscLogEventBegin(MATCOLORING_Comm, *etor, 0, 0, 0));
64 PetscCall(PetscSFComputeDegreeBegin(*etor, &rowdegrees));
65 PetscCall(PetscSFComputeDegreeEnd(*etor, &rowdegrees));
66 PetscCall(PetscLogEventEnd(MATCOLORING_Comm, *etor, 0, 0, 0));
67
68 PetscCall(PetscFree(rowdata));
69 PetscCall(PetscFree(rowleaf));
70 PetscCall(PetscFree(colleaf));
71 PetscFunctionReturn(PETSC_SUCCESS);
72 }
73