COIN-OR::LEMON - Graph Library

source: glpk-cmake/src/glpqmd.c @ 1:c445c931472f

Last change on this file since 1:c445c931472f was 1:c445c931472f, checked in by Alpar Juttner <alpar@…>, 10 years ago

Import glpk-4.45

  • Generated files and doc/notes are removed
File size: 17.9 KB
Line 
1/* glpqmd.c (quotient minimum degree algorithm) */
2
3/***********************************************************************
4*  This code is part of GLPK (GNU Linear Programming Kit).
5*
6*  THIS CODE IS THE RESULT OF TRANSLATION OF THE FORTRAN SUBROUTINES
7*  GENQMD, QMDRCH, QMDQT, QMDUPD, AND QMDMRG FROM THE BOOK:
8*
9*  ALAN GEORGE, JOSEPH W-H LIU. COMPUTER SOLUTION OF LARGE SPARSE
10*  POSITIVE DEFINITE SYSTEMS. PRENTICE-HALL, 1981.
11*
12*  THE TRANSLATION HAS BEEN DONE WITH THE PERMISSION OF THE AUTHORS
13*  OF THE ORIGINAL FORTRAN SUBROUTINES: ALAN GEORGE AND JOSEPH LIU,
14*  UNIVERSITY OF WATERLOO, WATERLOO, ONTARIO, CANADA.
15*
16*  The translation was made by Andrew Makhorin <mao@gnu.org>.
17*
18*  GLPK is free software: you can redistribute it and/or modify it
19*  under the terms of the GNU General Public License as published by
20*  the Free Software Foundation, either version 3 of the License, or
21*  (at your option) any later version.
22*
23*  GLPK is distributed in the hope that it will be useful, but WITHOUT
24*  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
25*  or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public
26*  License for more details.
27*
28*  You should have received a copy of the GNU General Public License
29*  along with GLPK. If not, see <http://www.gnu.org/licenses/>.
30***********************************************************************/
31
32#include "glpqmd.h"
33
34/***********************************************************************
35*  NAME
36*
37*  genqmd - GENeral Quotient Minimum Degree algorithm
38*
39*  SYNOPSIS
40*
41*  #include "glpqmd.h"
42*  void genqmd(int *neqns, int xadj[], int adjncy[], int perm[],
43*     int invp[], int deg[], int marker[], int rchset[], int nbrhd[],
44*     int qsize[], int qlink[], int *nofsub);
45*
46*  PURPOSE
47*
48*  This routine implements the minimum degree algorithm. It makes use
49*  of the implicit representation of the elimination graph by quotient
50*  graphs, and the notion of indistinguishable nodes.
51*
52*  CAUTION
53*
54*  The adjancy vector adjncy will be destroyed.
55*
56*  INPUT PARAMETERS
57*
58*  neqns  - number of equations;
59*  (xadj, adjncy) -
60*           the adjancy structure.
61*
62*  OUTPUT PARAMETERS
63*
64*  perm   - the minimum degree ordering;
65*  invp   - the inverse of perm.
66*
67*  WORKING PARAMETERS
68*
69*  deg    - the degree vector. deg[i] is negative means node i has been
70*           numbered;
71*  marker - a marker vector, where marker[i] is negative means node i
72*           has been merged with another nodeand thus can be ignored;
73*  rchset - vector used for the reachable set;
74*  nbrhd  - vector used for neighborhood set;
75*  qsize  - vector used to store the size of indistinguishable
76*           supernodes;
77*  qlink  - vector used to store indistinguishable nodes, i, qlink[i],
78*           qlink[qlink[i]], ... are the members of the supernode
79*           represented by i.
80*
81*  PROGRAM SUBROUTINES
82*
83*  qmdrch, qmdqt, qmdupd.
84***********************************************************************/
85
86void genqmd(int *_neqns, int xadj[], int adjncy[], int perm[],
87      int invp[], int deg[], int marker[], int rchset[], int nbrhd[],
88      int qsize[], int qlink[], int *_nofsub)
89{     int inode, ip, irch, j, mindeg, ndeg, nhdsze, node, np, num,
90         nump1, nxnode, rchsze, search, thresh;
91#     define neqns  (*_neqns)
92#     define nofsub (*_nofsub)
93      /* Initialize degree vector and other working variables. */
94      mindeg = neqns;
95      nofsub = 0;
96      for (node = 1; node <= neqns; node++)
97      {  perm[node] = node;
98         invp[node] = node;
99         marker[node] = 0;
100         qsize[node] = 1;
101         qlink[node] = 0;
102         ndeg = xadj[node+1] - xadj[node];
103         deg[node] = ndeg;
104         if (ndeg < mindeg) mindeg = ndeg;
105      }
106      num = 0;
107      /* Perform threshold search to get a node of min degree.
108         Variable search point to where search should start. */
109s200: search = 1;
110      thresh = mindeg;
111      mindeg = neqns;
112s300: nump1 = num + 1;
113      if (nump1 > search) search = nump1;
114      for (j = search; j <= neqns; j++)
115      {  node = perm[j];
116         if (marker[node] >= 0)
117         {  ndeg = deg[node];
118            if (ndeg <= thresh) goto s500;
119            if (ndeg < mindeg) mindeg = ndeg;
120         }
121      }
122      goto s200;
123      /* Node has minimum degree. Find its reachable sets by calling
124         qmdrch. */
125s500: search = j;
126      nofsub += deg[node];
127      marker[node] = 1;
128      qmdrch(&node, xadj, adjncy, deg, marker, &rchsze, rchset, &nhdsze,
129         nbrhd);
130      /* Eliminate all nodes indistinguishable from node. They are given
131         by node, qlink[node], ... . */
132      nxnode = node;
133s600: num++;
134      np = invp[nxnode];
135      ip = perm[num];
136      perm[np] = ip;
137      invp[ip] = np;
138      perm[num] = nxnode;
139      invp[nxnode] = num;
140      deg[nxnode] = -1;
141      nxnode = qlink[nxnode];
142      if (nxnode > 0) goto s600;
143      if (rchsze > 0)
144      {  /* Update the degrees of the nodes in the reachable set and
145            identify indistinguishable nodes. */
146         qmdupd(xadj, adjncy, &rchsze, rchset, deg, qsize, qlink,
147            marker, &rchset[rchsze+1], &nbrhd[nhdsze+1]);
148         /* Reset marker value of nodes in reach set. Update threshold
149            value for cyclic search. Also call qmdqt to form new
150            quotient graph. */
151         marker[node] = 0;
152         for (irch = 1; irch <= rchsze; irch++)
153         {  inode = rchset[irch];
154            if (marker[inode] >= 0)
155            {  marker[inode] = 0;
156               ndeg = deg[inode];
157               if (ndeg < mindeg) mindeg = ndeg;
158               if (ndeg <= thresh)
159               {  mindeg = thresh;
160                  thresh = ndeg;
161                  search = invp[inode];
162               }
163            }
164         }
165         if (nhdsze > 0)
166            qmdqt(&node, xadj, adjncy, marker, &rchsze, rchset, nbrhd);
167      }
168      if (num < neqns) goto s300;
169      return;
170#     undef neqns
171#     undef nofsub
172}
173
174/***********************************************************************
175*  NAME
176*
177*  qmdrch - Quotient MD ReaCHable set
178*
179*  SYNOPSIS
180*
181*  #include "glpqmd.h"
182*  void qmdrch(int *root, int xadj[], int adjncy[], int deg[],
183*     int marker[], int *rchsze, int rchset[], int *nhdsze,
184*     int nbrhd[]);
185*
186*  PURPOSE
187*
188*  This subroutine determines the reachable set of a node through a
189*  given subset. The adjancy structure is assumed to be stored in a
190*  quotient graph format.
191*
192*  INPUT PARAMETERS
193*
194*  root   - the given node not in the subset;
195*  (xadj, adjncy) -
196*           the adjancy structure pair;
197*  deg    - the degree vector. deg[i] < 0 means the node belongs to the
198*           given subset.
199*
200*  OUTPUT PARAMETERS
201*
202*  (rchsze, rchset) -
203*           the reachable set;
204*  (nhdsze, nbrhd) -
205*           the neighborhood set.
206*
207*  UPDATED PARAMETERS
208*
209*  marker - the marker vector for reach and nbrhd sets. > 0 means the
210*           node is in reach set. < 0 means the node has been merged
211*           with others in the quotient or it is in nbrhd set.
212***********************************************************************/
213
214void qmdrch(int *_root, int xadj[], int adjncy[], int deg[],
215      int marker[], int *_rchsze, int rchset[], int *_nhdsze,
216      int nbrhd[])
217{     int i, istop, istrt, j, jstop, jstrt, nabor, node;
218#     define root   (*_root)
219#     define rchsze (*_rchsze)
220#     define nhdsze (*_nhdsze)
221      /* Loop through the neighbors of root in the quotient graph. */
222      nhdsze = 0;
223      rchsze = 0;
224      istrt = xadj[root];
225      istop = xadj[root+1] - 1;
226      if (istop < istrt) return;
227      for (i = istrt; i <= istop; i++)
228      {  nabor = adjncy[i];
229         if (nabor == 0) return;
230         if (marker[nabor] == 0)
231         {  if (deg[nabor] >= 0)
232            {  /* Include nabor into the reachable set. */
233               rchsze++;
234               rchset[rchsze] = nabor;
235               marker[nabor] = 1;
236               goto s600;
237            }
238            /* nabor has been eliminated. Find nodes reachable from
239               it. */
240            marker[nabor] = -1;
241            nhdsze++;
242            nbrhd[nhdsze] = nabor;
243s300:       jstrt = xadj[nabor];
244            jstop = xadj[nabor+1] - 1;
245            for (j = jstrt; j <= jstop; j++)
246            {  node = adjncy[j];
247               nabor = - node;
248               if (node < 0) goto s300;
249               if (node == 0) goto s600;
250               if (marker[node] == 0)
251               {  rchsze++;
252                  rchset[rchsze] = node;
253                  marker[node] = 1;
254               }
255            }
256         }
257s600:    ;
258      }
259      return;
260#     undef root
261#     undef rchsze
262#     undef nhdsze
263}
264
265/***********************************************************************
266*  NAME
267*
268*  qmdqt - Quotient MD Quotient graph Transformation
269*
270*  SYNOPSIS
271*
272*  #include "glpqmd.h"
273*  void qmdqt(int *root, int xadj[], int adjncy[], int marker[],
274*     int *rchsze, int rchset[], int nbrhd[]);
275*
276*  PURPOSE
277*
278*  This subroutine performs the quotient graph transformation after a
279*  node has been eliminated.
280*
281*  INPUT PARAMETERS
282*
283*  root   - the node just eliminated. It becomes the representative of
284*           the new supernode;
285*  (xadj, adjncy) -
286*           the adjancy structure;
287*  (rchsze, rchset) -
288*           the reachable set of root in the old quotient graph;
289*  nbrhd  - the neighborhood set which will be merged with root to form
290*           the new supernode;
291*  marker - the marker vector.
292*
293*  UPDATED PARAMETERS
294*
295*  adjncy - becomes the adjncy of the quotient graph.
296***********************************************************************/
297
298void qmdqt(int *_root, int xadj[], int adjncy[], int marker[],
299      int *_rchsze, int rchset[], int nbrhd[])
300{     int inhd, irch, j, jstop, jstrt, link, nabor, node;
301#     define root   (*_root)
302#     define rchsze (*_rchsze)
303      irch = 0;
304      inhd = 0;
305      node = root;
306s100: jstrt = xadj[node];
307      jstop = xadj[node+1] - 2;
308      if (jstop >= jstrt)
309      {  /* Place reach nodes into the adjacent list of node. */
310         for (j = jstrt; j <= jstop; j++)
311         {  irch++;
312            adjncy[j] = rchset[irch];
313            if (irch >= rchsze) goto s400;
314         }
315      }
316      /* Link to other space provided by the nbrhd set. */
317      link = adjncy[jstop+1];
318      node = - link;
319      if (link >= 0)
320      {  inhd++;
321         node = nbrhd[inhd];
322         adjncy[jstop+1] = - node;
323      }
324      goto s100;
325      /* All reachable nodes have been saved. End the adjacent list.
326         Add root to the neighborhood list of each node in the reach
327         set. */
328s400: adjncy[j+1] = 0;
329      for (irch = 1; irch <= rchsze; irch++)
330      {  node = rchset[irch];
331         if (marker[node] >= 0)
332         {  jstrt = xadj[node];
333            jstop = xadj[node+1] - 1;
334            for (j = jstrt; j <= jstop; j++)
335            {  nabor = adjncy[j];
336               if (marker[nabor] < 0)
337               {  adjncy[j] = root;
338                  goto s600;
339               }
340            }
341         }
342s600:    ;
343      }
344      return;
345#     undef root
346#     undef rchsze
347}
348
349/***********************************************************************
350*  NAME
351*
352*  qmdupd - Quotient MD UPDate
353*
354*  SYNOPSIS
355*
356*  #include "glpqmd.h"
357*  void qmdupd(int xadj[], int adjncy[], int *nlist, int list[],
358*     int deg[], int qsize[], int qlink[], int marker[], int rchset[],
359*     int nbrhd[]);
360*
361*  PURPOSE
362*
363*  This routine performs degree update for a set of nodes in the minimum
364*  degree algorithm.
365*
366*  INPUT PARAMETERS
367*
368*  (xadj, adjncy) -
369*           the adjancy structure;
370*  (nlist, list) -
371*           the list of nodes whose degree has to be updated.
372*
373*  UPDATED PARAMETERS
374*
375*  deg    - the degree vector;
376*  qsize  - size of indistinguishable supernodes;
377*  qlink  - linked list for indistinguishable nodes;
378*  marker - used to mark those nodes in reach/nbrhd sets.
379*
380*  WORKING PARAMETERS
381*
382*  rchset - the reachable set;
383*  nbrhd  - the neighborhood set.
384*
385*  PROGRAM SUBROUTINES
386*
387*  qmdmrg.
388***********************************************************************/
389
390void qmdupd(int xadj[], int adjncy[], int *_nlist, int list[],
391      int deg[], int qsize[], int qlink[], int marker[], int rchset[],
392      int nbrhd[])
393{     int deg0, deg1, il, inhd, inode, irch, j, jstop, jstrt, mark,
394         nabor, nhdsze, node, rchsze;
395#     define nlist  (*_nlist)
396      /* Find all eliminated supernodes that are adjacent to some nodes
397         in the given list. Put them into (nhdsze, nbrhd). deg0 contains
398         the number of nodes in the list. */
399      if (nlist <= 0) return;
400      deg0 = 0;
401      nhdsze = 0;
402      for (il = 1; il <= nlist; il++)
403      {  node = list[il];
404         deg0 += qsize[node];
405         jstrt = xadj[node];
406         jstop = xadj[node+1] - 1;
407         for (j = jstrt; j <= jstop; j++)
408         {  nabor = adjncy[j];
409            if (marker[nabor] == 0 && deg[nabor] < 0)
410            {  marker[nabor] = -1;
411               nhdsze++;
412               nbrhd[nhdsze] = nabor;
413            }
414         }
415      }
416      /* Merge indistinguishable nodes in the list by calling the
417         subroutine qmdmrg. */
418      if (nhdsze > 0)
419         qmdmrg(xadj, adjncy, deg, qsize, qlink, marker, &deg0, &nhdsze,
420            nbrhd, rchset, &nbrhd[nhdsze+1]);
421      /* Find the new degrees of the nodes that have not been merged. */
422      for (il = 1; il <= nlist; il++)
423      {  node = list[il];
424         mark = marker[node];
425         if (mark == 0 || mark == 1)
426         {  marker[node] = 2;
427            qmdrch(&node, xadj, adjncy, deg, marker, &rchsze, rchset,
428               &nhdsze, nbrhd);
429            deg1 = deg0;
430            if (rchsze > 0)
431            {  for (irch = 1; irch <= rchsze; irch++)
432               {  inode = rchset[irch];
433                  deg1 += qsize[inode];
434                  marker[inode] = 0;
435               }
436            }
437            deg[node] = deg1 - 1;
438            if (nhdsze > 0)
439            {  for (inhd = 1; inhd <= nhdsze; inhd++)
440               {  inode = nbrhd[inhd];
441                  marker[inode] = 0;
442               }
443            }
444         }
445      }
446      return;
447#     undef nlist
448}
449
450/***********************************************************************
451*  NAME
452*
453*  qmdmrg - Quotient MD MeRGe
454*
455*  SYNOPSIS
456*
457*  #include "qmdmrg.h"
458*  void qmdmrg(int xadj[], int adjncy[], int deg[], int qsize[],
459*     int qlink[], int marker[], int *deg0, int *nhdsze, int nbrhd[],
460*     int rchset[], int ovrlp[]);
461*
462*  PURPOSE
463*
464*  This routine merges indistinguishable nodes in the minimum degree
465*  ordering algorithm. It also computes the new degrees of these new
466*  supernodes.
467*
468*  INPUT PARAMETERS
469*
470*  (xadj, adjncy) -
471*           the adjancy structure;
472*  deg0   - the number of nodes in the given set;
473*  (nhdsze, nbrhd) -
474*           the set of eliminated supernodes adjacent to some nodes in
475*           the set.
476*
477*  UPDATED PARAMETERS
478*
479*  deg    - the degree vector;
480*  qsize  - size of indistinguishable nodes;
481*  qlink  - linked list for indistinguishable nodes;
482*  marker - the given set is given by those nodes with marker value set
483*           to 1. Those nodes with degree updated will have marker value
484*           set to 2.
485*
486*  WORKING PARAMETERS
487*
488*  rchset - the reachable set;
489*  ovrlp  - temp vector to store the intersection of two reachable sets.
490***********************************************************************/
491
492void qmdmrg(int xadj[], int adjncy[], int deg[], int qsize[],
493      int qlink[], int marker[], int *_deg0, int *_nhdsze, int nbrhd[],
494      int rchset[], int ovrlp[])
495{     int deg1, head, inhd, iov, irch, j, jstop, jstrt, link, lnode,
496         mark, mrgsze, nabor, node, novrlp, rchsze, root;
497#     define deg0   (*_deg0)
498#     define nhdsze (*_nhdsze)
499      /* Initialization. */
500      if (nhdsze <= 0) return;
501      for (inhd = 1; inhd <= nhdsze; inhd++)
502      {  root = nbrhd[inhd];
503         marker[root] = 0;
504      }
505      /* Loop through each eliminated supernode in the set
506         (nhdsze, nbrhd). */
507      for (inhd = 1; inhd <= nhdsze; inhd++)
508      {  root = nbrhd[inhd];
509         marker[root] = -1;
510         rchsze = 0;
511         novrlp = 0;
512         deg1 = 0;
513s200:    jstrt = xadj[root];
514         jstop = xadj[root+1] - 1;
515         /* Determine the reachable set and its intersection with the
516            input reachable set. */
517         for (j = jstrt; j <= jstop; j++)
518         {  nabor = adjncy[j];
519            root = - nabor;
520            if (nabor < 0) goto s200;
521            if (nabor == 0) break;
522            mark = marker[nabor];
523            if (mark == 0)
524            {  rchsze++;
525               rchset[rchsze] = nabor;
526               deg1 += qsize[nabor];
527               marker[nabor] = 1;
528            }
529            else if (mark == 1)
530            {  novrlp++;
531               ovrlp[novrlp] = nabor;
532               marker[nabor] = 2;
533            }
534         }
535         /* From the overlapped set, determine the nodes that can be
536            merged together. */
537         head = 0;
538         mrgsze = 0;
539         for (iov = 1; iov <= novrlp; iov++)
540         {  node = ovrlp[iov];
541            jstrt = xadj[node];
542            jstop = xadj[node+1] - 1;
543            for (j = jstrt; j <= jstop; j++)
544            {  nabor = adjncy[j];
545               if (marker[nabor] == 0)
546               {  marker[node] = 1;
547                  goto s1100;
548               }
549            }
550            /* Node belongs to the new merged supernode. Update the
551               vectors qlink and qsize. */
552            mrgsze += qsize[node];
553            marker[node] = -1;
554            lnode = node;
555s900:       link = qlink[lnode];
556            if (link > 0)
557            {  lnode = link;
558               goto s900;
559            }
560            qlink[lnode] = head;
561            head = node;
562s1100:      ;
563         }
564         if (head > 0)
565         {  qsize[head] = mrgsze;
566            deg[head] = deg0 + deg1 - 1;
567            marker[head] = 2;
568         }
569         /* Reset marker values. */
570         root = nbrhd[inhd];
571         marker[root] = 0;
572         if (rchsze > 0)
573         {  for (irch = 1; irch <= rchsze; irch++)
574            {  node = rchset[irch];
575               marker[node] = 0;
576            }
577         }
578      }
579      return;
580#     undef deg0
581#     undef nhdsze
582}
583
584/* eof */
Note: See TracBrowser for help on using the repository browser.