1 /* ido.f -- translated by f2c (version of 25 March 1992 12:58:56).*/
2
3 #include <../src/mat/graphops/color/impls/minpack/color.h>
4
5 static PetscInt c_n1 = -1;
6
MINPACKido(PetscInt * m,PetscInt * n,const PetscInt * indrow,const PetscInt * jpntr,const PetscInt * indcol,const PetscInt * ipntr,PetscInt * ndeg,PetscInt * list,PetscInt * maxclq,PetscInt * iwa1,PetscInt * iwa2,PetscInt * iwa3,PetscInt * iwa4)7 PetscErrorCode MINPACKido(PetscInt *m, PetscInt *n, const PetscInt *indrow, const PetscInt *jpntr, const PetscInt *indcol, const PetscInt *ipntr, PetscInt *ndeg, PetscInt *list, PetscInt *maxclq, PetscInt *iwa1, PetscInt *iwa2, PetscInt *iwa3, PetscInt *iwa4)
8 {
9 /* System generated locals */
10 PetscInt i__1, i__2, i__3, i__4;
11
12 /* Local variables */
13 PetscInt jcol = 0, ncomp = 0, ic, ip, jp, ir, maxinc, numinc, numord, maxlst, numwgt, numlst;
14
15 /* Given the sparsity pattern of an m by n matrix A, this */
16 /* subroutine determines an incidence-degree ordering of the */
17 /* columns of A. */
18 /* The incidence-degree ordering is defined for the loopless */
19 /* graph G with vertices a(j), j = 1,2,...,n where a(j) is the */
20 /* j-th column of A and with edge (a(i),a(j)) if and only if */
21 /* columns i and j have a non-zero in the same row position. */
22 /* The incidence-degree ordering is determined recursively by */
23 /* letting list(k), k = 1,...,n be a column with maximal */
24 /* incidence to the subgraph spanned by the ordered columns. */
25 /* Among all the columns of maximal incidence, ido chooses a */
26 /* column of maximal degree. */
27 /* The subroutine statement is */
28 /* subroutine ido(m,n,indrow,jpntr,indcol,ipntr,ndeg,list, */
29 /* maxclq,iwa1,iwa2,iwa3,iwa4) */
30 /* where */
31 /* m is a positive integer input variable set to the number */
32 /* of rows of A. */
33 /* n is a positive integer input variable set to the number */
34 /* of columns of A. */
35 /* indrow is an integer input array which contains the row */
36 /* indices for the non-zeroes in the matrix A. */
37 /* jpntr is an integer input array of length n + 1 which */
38 /* specifies the locations of the row indices in indrow. */
39 /* The row indices for column j are */
40 /* indrow(k), k = jpntr(j),...,jpntr(j+1)-1. */
41 /* Note that jpntr(n+1)-1 is then the number of non-zero */
42 /* elements of the matrix A. */
43 /* indcol is an integer input array which contains the */
44 /* column indices for the non-zeroes in the matrix A. */
45 /* ipntr is an integer input array of length m + 1 which */
46 /* specifies the locations of the column indices in indcol. */
47 /* The column indices for row i are */
48 /* indcol(k), k = ipntr(i),...,ipntr(i+1)-1. */
49 /* Note that ipntr(m+1)-1 is then the number of non-zero */
50 /* elements of the matrix A. */
51 /* ndeg is an integer input array of length n which specifies */
52 /* the degree sequence. The degree of the j-th column */
53 /* of A is ndeg(j). */
54 /* list is an integer output array of length n which specifies */
55 /* the incidence-degree ordering of the columns of A. The j-th */
56 /* column in this order is list(j). */
57 /* maxclq is an integer output variable set to the size */
58 /* of the largest clique found during the ordering. */
59 /* iwa1,iwa2,iwa3, and iwa4 are integer work arrays of length n. */
60 /* Subprograms called */
61 /* MINPACK-supplied ... numsrt */
62 /* FORTRAN-supplied ... max */
63 /* Argonne National Laboratory. MINPACK Project. August 1984. */
64 /* Thomas F. Coleman, Burton S. Garbow, Jorge J. More' */
65
66 /* Sort the degree sequence. */
67
68 PetscFunctionBegin;
69 /* Parameter adjustments */
70 --iwa4;
71 --iwa3;
72 --iwa2;
73 --list;
74 --ndeg;
75 --ipntr;
76 --indcol;
77 --jpntr;
78 --indrow;
79
80 /* Function Body */
81 i__1 = *n - 1;
82 PetscCall(MINPACKnumsrt(n, &i__1, &ndeg[1], &c_n1, &iwa4[1], &iwa2[1], &iwa3[1]));
83
84 /* Initialization block. */
85 /* Create a doubly-linked list to access the incidences of the */
86 /* columns. The pointers for the linked list are as follows. */
87 /* Each un-ordered column ic is in a list (the incidence list) */
88 /* of columns with the same incidence. */
89 /* iwa1(numinc) is the first column in the numinc list */
90 /* unless iwa1(numinc) = 0. In this case there are */
91 /* no columns in the numinc list. */
92 /* iwa2(ic) is the column before ic in the incidence list */
93 /* unless iwa2(ic) = 0. In this case ic is the first */
94 /* column in this incidence list. */
95 /* iwa3(ic) is the column after ic in the incidence list */
96 /* unless iwa3(ic) = 0. In this case ic is the last */
97 /* column in this incidence list. */
98 /* If ic is an un-ordered column, then list(ic) is the */
99 /* incidence of ic to the graph induced by the ordered */
100 /* columns. If jcol is an ordered column, then list(jcol) */
101 /* is the incidence-degree order of column jcol. */
102
103 maxinc = 0;
104 for (jp = *n; jp >= 1; --jp) {
105 ic = iwa4[jp];
106 iwa1[*n - jp] = 0;
107 iwa2[ic] = 0;
108 iwa3[ic] = iwa1[0];
109 if (iwa1[0] > 0) iwa2[iwa1[0]] = ic;
110 iwa1[0] = ic;
111 iwa4[jp] = 0;
112 list[jp] = 0;
113 }
114
115 /* Determine the maximal search length for the list */
116 /* of columns of maximal incidence. */
117
118 maxlst = 0;
119 i__1 = *m;
120 for (ir = 1; ir <= i__1; ++ir) {
121 /* Computing 2nd power */
122 i__2 = ipntr[ir + 1] - ipntr[ir];
123 maxlst += i__2 * i__2;
124 }
125 maxlst /= *n;
126 *maxclq = 0;
127 numord = 1;
128
129 /* Beginning of iteration loop. */
130
131 L30:
132
133 /* Choose a column jcol of maximal degree among the */
134 /* columns of maximal incidence maxinc. */
135
136 L40:
137 jp = iwa1[maxinc];
138 if (jp > 0) goto L50;
139 --maxinc;
140 goto L40;
141 L50:
142 numwgt = -1;
143 i__1 = maxlst;
144 for (numlst = 1; numlst <= i__1; ++numlst) {
145 if (ndeg[jp] > numwgt) {
146 numwgt = ndeg[jp];
147 jcol = jp;
148 }
149 jp = iwa3[jp];
150 if (jp <= 0) goto L70;
151 }
152 L70:
153 list[jcol] = numord;
154
155 /* Update the size of the largest clique */
156 /* found during the ordering. */
157
158 if (!maxinc) ncomp = 0;
159 ++ncomp;
160 if (maxinc + 1 == ncomp) *maxclq = PetscMax(*maxclq, ncomp);
161
162 /* Termination test. */
163
164 ++numord;
165 if (numord > *n) goto L100;
166
167 /* Delete column jcol from the maxinc list. */
168
169 if (!iwa2[jcol]) iwa1[maxinc] = iwa3[jcol];
170 else iwa3[iwa2[jcol]] = iwa3[jcol];
171
172 if (iwa3[jcol] > 0) iwa2[iwa3[jcol]] = iwa2[jcol];
173
174 /* Find all columns adjacent to column jcol. */
175
176 iwa4[jcol] = *n;
177
178 /* Determine all positions (ir,jcol) which correspond */
179 /* to non-zeroes in the matrix. */
180
181 i__1 = jpntr[jcol + 1] - 1;
182 for (jp = jpntr[jcol]; jp <= i__1; ++jp) {
183 ir = indrow[jp];
184
185 /* For each row ir, determine all positions (ir,ic) */
186 /* which correspond to non-zeroes in the matrix. */
187
188 i__2 = ipntr[ir + 1] - 1;
189 for (ip = ipntr[ir]; ip <= i__2; ++ip) {
190 ic = indcol[ip];
191
192 /* Array iwa4 marks columns which are adjacent to */
193 /* column jcol. */
194
195 if (iwa4[ic] < numord) {
196 iwa4[ic] = numord;
197
198 /* Update the pointers to the current incidence lists. */
199
200 numinc = list[ic];
201 ++list[ic];
202 /* Computing MAX */
203 i__3 = maxinc;
204 i__4 = list[ic];
205 maxinc = PetscMax(i__3, i__4);
206
207 /* Delete column ic from the numinc list. */
208
209 if (!iwa2[ic]) iwa1[numinc] = iwa3[ic];
210 else iwa3[iwa2[ic]] = iwa3[ic];
211
212 if (iwa3[ic] > 0) iwa2[iwa3[ic]] = iwa2[ic];
213
214 /* Add column ic to the numinc+1 list. */
215
216 iwa2[ic] = 0;
217 iwa3[ic] = iwa1[numinc + 1];
218 if (iwa1[numinc + 1] > 0) iwa2[iwa1[numinc + 1]] = ic;
219 iwa1[numinc + 1] = ic;
220 }
221 }
222 }
223
224 /* End of iteration loop. */
225
226 goto L30;
227 L100:
228
229 /* Invert the array list. */
230
231 i__1 = *n;
232 for (jcol = 1; jcol <= i__1; ++jcol) iwa2[list[jcol]] = jcol;
233 i__1 = *n;
234 for (jp = 1; jp <= i__1; ++jp) list[jp] = iwa2[jp];
235 PetscFunctionReturn(PETSC_SUCCESS);
236 }
237