symbfct.c

斯坦福大学Grant和Boyd教授等开发的凸优化matlab工具箱

    /* Local variables */
    static integer temp, xsup, last1, last2, j, k, lflag, pleaf, hinbr, jstop,
	     jstrt, ifdesc, oldnbr, parent, lownbr, lca;


/*       ----------- */
/*       PARAMETERS. */
/*       ----------- */

/*       ---------------- */
/*       LOCAL VARIABLES. */
/*       ---------------- */

/* ***********************************************************************
 */

/*       -------------------------------------------------- */
/*       COMPUTE LEVEL(*), FDESC(*), NCHILD(*). */
/*       INITIALIZE ROWCNT(*), COLCNT(*), */
/*                  SET(*), PRVLF(*), WEIGHT(*), PRVNBR(*). */
/*       -------------------------------------------------- */
    /* Parameter adjustments */
    --prvnbr;
    --prvlf;
    --set;
    --colcnt;
    --rowcnt;
    --etpar;
    --invp;
    --perm;
    --adjncy;
    --xadj;

    /* Function Body */
    level[0] = 0;
    for (k = *neqns; k >= 1; --k) {
	rowcnt[k] = 1;
	colcnt[k] = 0;
	set[k] = k;
	prvlf[k] = 0;
	level[k] = level[etpar[k]] + 1;
	weight[k] = 1;
	fdesc[k] = k;
	nchild[k] = 0;
	prvnbr[k] = 0;
/* L100: */
    }
    nchild[0] = 0;
    fdesc[0] = 0;
    i__1 = *neqns;
    for (k = 1; k <= i__1; ++k) {
	parent = etpar[k];
	weight[parent] = 0;
	++nchild[parent];
	ifdesc = fdesc[k];
	if (ifdesc < fdesc[parent]) {
	    fdesc[parent] = ifdesc;
	}
/* L200: */
    }
/*       ------------------------------------ */
/*       FOR EACH ``LOW NEIGHBOR'' LOWNBR ... */
/*       ------------------------------------ */
    i__1 = *neqns;
    for (lownbr = 1; lownbr <= i__1; ++lownbr) {
	lflag = 0;
	ifdesc = fdesc[lownbr];
	oldnbr = perm[lownbr];
	jstrt = xadj[oldnbr];
	jstop = xadj[oldnbr + 1] - 1;
/*           ----------------------------------------------- */
/*           FOR EACH ``HIGH NEIGHBOR'', HINBR OF LOWNBR ... */
/*           ----------------------------------------------- */
	i__2 = jstop;
	for (j = jstrt; j <= i__2; ++j) {
	    hinbr = invp[adjncy[j]];
	    if (hinbr > lownbr) {
		if (ifdesc > prvnbr[hinbr]) {
/*                       ------------------------- */
/*                       INCREMENT WEIGHT(LOWNBR). */
/*                       ------------------------- */
		    ++weight[lownbr];
		    pleaf = prvlf[hinbr];
/*                       --------------------------------
--------- */
/*                       IF HINBR HAS NO PREVIOUS ``LOW NE
IGHBOR'' */
/*                       THEN ... */
/*                       --------------------------------
--------- */
		    if (pleaf == 0) {
/*                           ------------------------
----------------- */
/*                           ... ACCUMULATE LOWNBR-->H
INBR PATH LENGTH */
/*                               IN ROWCNT(HINBR). */
/*                           ------------------------
----------------- */
			rowcnt[hinbr] = rowcnt[hinbr] + level[lownbr] - level[
				hinbr];
		    } else {
/*                           ------------------------
----------------- */
/*                           ... OTHERWISE, LCA <-- FI
ND(PLEAF), WHICH */
/*                               IS THE LEAST COMMON A
NCESTOR OF PLEAF */
/*                               AND LOWNBR. */
/*                               (PATH HALVING.) */
/*                           ------------------------
----------------- */
			last1 = pleaf;
			last2 = set[last1];
			lca = set[last2];
L300:
			if (lca != last2) {
			    set[last1] = lca;
			    last1 = lca;
			    last2 = set[last1];
			    lca = set[last2];
			    goto L300;
			}
/*                           ------------------------
------------- */
/*                           ACCUMULATE PLEAF-->LCA PA
TH LENGTH IN */
/*                           ROWCNT(HINBR). */
/*                           DECREMENT WEIGHT(LCA). */
/*                           ------------------------
------------- */
			rowcnt[hinbr] = rowcnt[hinbr] + level[lownbr] - level[
				lca];
			--weight[lca];
		    }
/*                       --------------------------------
-------------- */
/*                       LOWNBR NOW BECOMES ``PREVIOUS LEA
F'' OF HINBR. */
/*                       --------------------------------
-------------- */
		    prvlf[hinbr] = lownbr;
		    lflag = 1;
		}
/*                   ----------------------------------------
---------- */
/*                   LOWNBR NOW BECOMES ``PREVIOUS NEIGHBOR'' 
OF HINBR. */
/*                   ----------------------------------------
---------- */
		prvnbr[hinbr] = lownbr;
	    }
/* L500: */
	}
/*           ---------------------------------------------------- */
/*           DECREMENT WEIGHT ( PARENT(LOWNBR) ). */
/*           SET ( P(LOWNBR) ) <-- SET ( P(LOWNBR) ) + SET(XSUP). */
/*           ---------------------------------------------------- */
	parent = etpar[lownbr];
	--weight[parent];
	if (lflag == 1 || nchild[lownbr] >= 2) {
	    xsup = lownbr;
	}
	set[xsup] = parent;
/* L600: */
    }
/*       --------------------------------------------------------- */
/*       USE WEIGHTS TO COMPUTE COLUMN (AND TOTAL) NONZERO COUNTS. */
/*       --------------------------------------------------------- */
    *nlnz = 0;
    i__1 = *neqns;
    for (k = 1; k <= i__1; ++k) {
	temp = colcnt[k] + weight[k];
	colcnt[k] = temp;
	*nlnz += temp;
	parent = etpar[k];
	if (parent != 0) {
	    colcnt[parent] += temp;
	}
/* L700: */
    }

    return 0;
} /* fcnthn_ */

/* *********************************************************************** */
/* *********************************************************************** */

/*   Version:        0.3 */
/*   Last modified:  December 27, 1994 */
/*   Authors:        Esmond G. Ng and Barry W. Peyton */

/*   Mathematical Sciences Section, Oak Ridge National Laboratory */

/* *********************************************************************** */
/* *********************************************************************** */
/* ****************    FSUP1 ..... FIND SUPERNODES #1    ***************** */
/* *********************************************************************** */
/* *********************************************************************** */

/*   PURPOSE: */
/*       THIS SUBROUTINE IS THE FIRST OF TWO ROUTINES FOR FINDING A */
/*       MAXIMAL SUPERNODE PARTITION.  IT RETURNS ONLY THE NUMBER OF */
/*       SUPERNODES NSUPER AND THE SUPERNODE MEMBERSHIP VECTOR SNODE(*), */
/*       WHICH IS OF LENGTH NEQNS.  THE VECTORS OF LENGTH NSUPER ARE */
/*       COMPUTED SUBSEQUENTLY BY THE COMPANION ROUTINE FSUP2. */

/*   METHOD AND ASSUMPTIONS: */
/*       THIS ROUTINE USES THE ELIMINATION TREE AND THE FACTOR COLUMN */
/*       COUNTS TO COMPUTE THE SUPERNODE PARTITION; IT ALSO ASSUMES A */
/*       POSTORDERING OF THE ELIMINATION TREE. */

/*   INPUT PARAMETERS: */
/*       (I) NEQNS       -   NUMBER OF EQUATIONS. */
/*       (I) ETPAR(*)    -   ARRAY OF LENGTH NEQNS, CONTAINING THE */
/*                           ELIMINATION TREE OF THE POSTORDERED MATRIX. */
/*       (I) COLCNT(*)   -   ARRAY OF LENGTH NEQNS, CONTAINING THE */
/*                           FACTOR COLUMN COUNTS: I.E., THE NUMBER OF */
/*                           NONZERO ENTRIES IN EACH COLUMN OF L */
/*                           (INCLUDING THE DIAGONAL ENTRY). */

/*   OUTPUT PARAMETERS: */
/*       (I) NOFSUB      -   NUMBER OF SUBSCRIPTS. */
/*       (I) NSUPER      -   NUMBER OF SUPERNODES (<= NEQNS). */
/*       (I) SNODE(*)    -   ARRAY OF LENGTH NEQNS FOR RECORDING */
/*                           SUPERNODE MEMBERSHIP. */

/*   FIRST CREATED ON    JANUARY 18, 1992. */
/*   LAST UPDATED ON     NOVEMBER 11, 1994. */

/* *********************************************************************** */

/* Subroutine */ int fsup1_(neqns, etpar, colcnt, nofsub, nsuper, snode)
integer *neqns, *etpar, *colcnt, *nofsub, *nsuper, *snode;
{
    /* System generated locals */
    integer i__1;

    /* Local variables */
    static integer kcol;


/* ***********************************************************************
 */

/*       ----------- */
/*       PARAMETERS. */
/*       ----------- */

/*       ---------------- */
/*       LOCAL VARIABLES. */
/*       ---------------- */

/* ***********************************************************************
 */

/*       -------------------------------------------- */
/*       COMPUTE THE FUNDAMENTAL SUPERNODE PARTITION. */
/*       -------------------------------------------- */
    /* Parameter adjustments */
    --snode;
    --colcnt;
    --etpar;

    /* Function Body */
    *nsuper = 1;
    snode[1] = 1;
    *nofsub = colcnt[1];
    i__1 = *neqns;
    for (kcol = 2; kcol <= i__1; ++kcol) {
	if (etpar[kcol - 1] == kcol) {
	    if (colcnt[kcol - 1] == colcnt[kcol] + 1) {
		snode[kcol] = *nsuper;
		goto L300;
	    }
	}
	++(*nsuper);
	snode[kcol] = *nsuper;
	*nofsub += colcnt[kcol];
L300:
	;
    }

    return 0;
} /* fsup1_ */

/* *********************************************************************** */
/* *********************************************************************** */

/*   Version:        0.3 */
/*   Last modified:  December 27, 1994 */
/*   Authors:        Esmond G. Ng and Barry W. Peyton */

/*   Mathematical Sciences Section, Oak Ridge National Laboratory */

/* *********************************************************************** */
/* *********************************************************************** */
/* ****************    FSUP2  ..... FIND SUPERNODES #2   ***************** */
/* *********************************************************************** */
/* *********************************************************************** */

/*   PURPOSE: */
/*       THIS SUBROUTINE IS THE SECOND OF TWO ROUTINES FOR FINDING A */
/*       MAXIMAL SUPERNODE PARTITION.  IT'S SOLE PURPOSE IS TO */
/*       CONSTRUCT THE NEEDED VECTOR OF LENGTH NSUPER: XSUPER(*).  THE */
/*       FIRST ROUTINE FSUP1 COMPUTES THE NUMBER OF SUPERNODES AND THE */
/*       SUPERNODE MEMBERSHIP VECTOR SNODE(*), WHICH IS OF LENGTH NEQNS. */


/*   ASSUMPTIONS: */
/*       THIS ROUTINE ASSUMES A POSTORDERING OF THE ELIMINATION TREE.  IT */
/*       ALSO ASSUMES THAT THE OUTPUT FROM FSUP1 IS AVAILABLE. */

/*   INPUT PARAMETERS: */
/*       (I) NEQNS       -   NUMBER OF EQUATIONS. */
/*       (I) NSUPER      -   NUMBER OF SUPERNODES (<= NEQNS). */
/*       (I) ETPAR(*)    -   ARRAY OF LENGTH NEQNS, CONTAINING THE */
/*                           ELIMINATION TREE OF THE POSTORDERED MATRIX. */
/*       (I) SNODE(*)    -   ARRAY OF LENGTH NEQNS FOR RECORDING */
/*                           SUPERNODE MEMBERSHIP. */

/*   OUTPUT PARAMETERS: */
/*       (I) XSUPER(*)   -   ARRAY OF LENGTH NSUPER+1, CONTAINING THE */
/*                           SUPERNODE PARTITIONING. */

/*   FIRST CREATED ON    JANUARY 18, 1992. */
/*   LAST UPDATED ON     NOVEMEBER 22, 1994. */

/* *********************************************************************** */

/* Subroutine */ int fsup2_(neqns, nsuper, etpar, snode, xsuper)
integer *neqns, *nsuper, *etpar, *snode, *xsuper;
{
    static integer kcol, ksup, lstsup;


/* ***********************************************************************
 */

/*       ----------- */
/*       PARAMETERS. */
/*       ----------- */

/*       ---------------- */
/*       LOCAL VARIABLES. */
/*       ---------------- */

/* ***********************************************************************
 */

/*       ------------------------------------------------- */
/*       COMPUTE THE SUPERNODE PARTITION VECTOR XSUPER(*). */
/*       ------------------------------------------------- */
    /* Parameter adjustments */
    --xsuper;
    --snode;
    --etpar;

    /* Function Body */
    lstsup = *nsuper + 1;
    for (kcol = *neqns; kcol >= 1; --kcol) {
	ksup = snode[kcol];
	if (ksup != lstsup) {
	    xsuper[lstsup] = kcol + 1;
	}
	lstsup = ksup;
/* L100: */
    }
    xsuper[1] = 1;

    return 0;
} /* fsup2_ */

/* *********************************************************************** */
/* *********************************************************************** */

/*   Version:        0.3 */
/*   Last modified:  February 13, 1995 */
/*   Authors:        Esmond G. Ng and Barry W. Peyton */

/*   Mathematical Sciences Section, Oak Ridge National Laboratory */

/* *********************************************************************** */
/* *********************************************************************** */
/* *************     SYMFCT ..... SYMBOLIC FACTORIZATION    ************** */
/* *********************************************************************** */
/* *********************************************************************** */