1 /* degr.f -- translated by f2c (version of 25 March 1992 12:58:56). */
2
3 #include <../src/mat/graphops/color/impls/minpack/color.h>
4
MINPACKdegr(PetscInt * n,const PetscInt * indrow,const PetscInt * jpntr,const PetscInt * indcol,const PetscInt * ipntr,PetscInt * ndeg,PetscInt * iwa)5 PetscErrorCode MINPACKdegr(PetscInt *n, const PetscInt *indrow, const PetscInt *jpntr, const PetscInt *indcol, const PetscInt *ipntr, PetscInt *ndeg, PetscInt *iwa)
6 {
7 /* System generated locals */
8 PetscInt i__1, i__2, i__3;
9
10 /* Local variables */
11 PetscInt jcol, ic, ip, jp, ir;
12
13 /* subroutine degr */
14 /* Given the sparsity pattern of an m by n matrix A, */
15 /* this subroutine determines the degree sequence for */
16 /* the intersection graph of the columns of A. */
17 /* In graph-theory terminology, the intersection graph of */
18 /* the columns of A is the loopless graph G with vertices */
19 /* a(j), j = 1,2,...,n where a(j) is the j-th column of A */
20 /* and with edge (a(i),a(j)) if and only if columns i and j */
21 /* have a non-zero in the same row position. */
22 /* Note that the value of m is not needed by degr and is */
23 /* therefore not present in the subroutine statement. */
24 /* The subroutine statement is */
25 /* subroutine degr(n,indrow,jpntr,indcol,ipntr,ndeg,iwa) */
26 /* where */
27 /* n is a positive integer input variable set to the number */
28 /* of columns of A. */
29 /* indrow is an integer input array which contains the row */
30 /* indices for the non-zeroes in the matrix A. */
31 /* jpntr is an integer input array of length n + 1 which */
32 /* specifies the locations of the row indices in indrow. */
33 /* The row indices for column j are */
34 /* indrow(k), k = jpntr(j),...,jpntr(j+1)-1. */
35 /* Note that jpntr(n+1)-1 is then the number of non-zero */
36 /* elements of the matrix A. */
37 /* indcol is an integer input array which contains the */
38 /* column indices for the non-zeroes in the matrix A. */
39 /* ipntr is an integer input array of length m + 1 which */
40 /* specifies the locations of the column indices in indcol. */
41 /* The column indices for row i are */
42 /* indcol(k), k = ipntr(i),...,ipntr(i+1)-1. */
43 /* Note that ipntr(m+1)-1 is then the number of non-zero */
44 /* elements of the matrix A. */
45 /* ndeg is an integer output array of length n which */
46 /* specifies the degree sequence. The degree of the */
47 /* j-th column of A is ndeg(j). */
48 /* iwa is an integer work array of length n. */
49 /* Argonne National Laboratory. MINPACK Project. July 1983. */
50 /* Thomas F. Coleman, Burton S. Garbow, Jorge J. More' */
51
52 PetscFunctionBegin;
53 /* Parameter adjustments */
54 --iwa;
55 --ndeg;
56 --ipntr;
57 --indcol;
58 --jpntr;
59 --indrow;
60
61 /* Function Body */
62 i__1 = *n;
63 for (jp = 1; jp <= i__1; ++jp) {
64 ndeg[jp] = 0;
65 iwa[jp] = 0;
66 }
67
68 /* Compute the degree sequence by determining the contributions */
69 /* to the degrees from the current(jcol) column and further */
70 /* columns which have not yet been considered. */
71
72 i__1 = *n;
73 for (jcol = 2; jcol <= i__1; ++jcol) {
74 iwa[jcol] = *n;
75
76 /* Determine all positions (ir,jcol) which correspond */
77 /* to non-zeroes in the matrix. */
78
79 i__2 = jpntr[jcol + 1] - 1;
80 for (jp = jpntr[jcol]; jp <= i__2; ++jp) {
81 ir = indrow[jp];
82
83 /* For each row ir, determine all positions (ir,ic) */
84 /* which correspond to non-zeroes in the matrix. */
85
86 i__3 = ipntr[ir + 1] - 1;
87 for (ip = ipntr[ir]; ip <= i__3; ++ip) {
88 ic = indcol[ip];
89
90 /* Array iwa marks columns which have contributed to */
91 /* the degree count of column jcol. Update the degree */
92 /* counts of these columns as well as column jcol. */
93
94 if (iwa[ic] < jcol) {
95 iwa[ic] = jcol;
96 ++ndeg[ic];
97 ++ndeg[jcol];
98 }
99 }
100 }
101 }
102 PetscFunctionReturn(PETSC_SUCCESS);
103 }
104