kusano 7d535a
/*! @file ilu_heap_relax_snode.c
kusano 7d535a
 * \brief Identify the initial relaxed supernodes
kusano 7d535a
 *
kusano 7d535a
 * 
kusano 7d535a
 * -- SuperLU routine (version 4.0) --
kusano 7d535a
 * Lawrence Berkeley National Laboratory
kusano 7d535a
 * June 1, 2009
kusano 7d535a
 * 
kusano 7d535a
 */
kusano 7d535a
kusano 7d535a
#include "slu_ddefs.h"
kusano 7d535a
kusano 7d535a
/*! \brief
kusano 7d535a
 *
kusano 7d535a
 * 
kusano 7d535a
 * Purpose
kusano 7d535a
 * =======
kusano 7d535a
 *    ilu_heap_relax_snode() - Identify the initial relaxed supernodes,
kusano 7d535a
 *    assuming that the matrix has been reordered according to the postorder
kusano 7d535a
 *    of the etree.
kusano 7d535a
 * 
kusano 7d535a
 */
kusano 7d535a
kusano 7d535a
void
kusano 7d535a
ilu_heap_relax_snode (
kusano 7d535a
	     const     int n,
kusano 7d535a
	     int       *et,	      /* column elimination tree */
kusano 7d535a
	     const int relax_columns, /* max no of columns allowed in a
kusano 7d535a
					 relaxed snode */
kusano 7d535a
	     int       *descendants,  /* no of descendants of each node
kusano 7d535a
					 in the etree */
kusano 7d535a
	     int       *relax_end,    /* last column in a supernode
kusano 7d535a
				       * if j-th column starts a relaxed
kusano 7d535a
				       * supernode, relax_end[j] represents
kusano 7d535a
				       * the last column of this supernode */
kusano 7d535a
	     int       *relax_fsupc   /* first column in a supernode
kusano 7d535a
				       * relax_fsupc[j] represents the first
kusano 7d535a
				       * column of j-th supernode */
kusano 7d535a
	     )
kusano 7d535a
{
kusano 7d535a
    register int i, j, k, l, f, parent;
kusano 7d535a
    register int snode_start;	/* beginning of a snode */
kusano 7d535a
    int *et_save, *post, *inv_post, *iwork;
kusano 7d535a
    int nsuper_et = 0, nsuper_et_post = 0;
kusano 7d535a
kusano 7d535a
    /* The etree may not be postordered, but is heap ordered. */
kusano 7d535a
kusano 7d535a
    iwork = (int*) intMalloc(3*n+2);
kusano 7d535a
    if ( !iwork ) ABORT("SUPERLU_MALLOC fails for iwork[]");
kusano 7d535a
    inv_post = iwork + n+1;
kusano 7d535a
    et_save = inv_post + n+1;
kusano 7d535a
kusano 7d535a
    /* Post order etree */
kusano 7d535a
    post = (int *) TreePostorder(n, et);
kusano 7d535a
    for (i = 0; i < n+1; ++i) inv_post[post[i]] = i;
kusano 7d535a
kusano 7d535a
    /* Renumber etree in postorder */
kusano 7d535a
    for (i = 0; i < n; ++i) {
kusano 7d535a
	iwork[post[i]] = post[et[i]];
kusano 7d535a
	et_save[i] = et[i]; /* Save the original etree */
kusano 7d535a
    }
kusano 7d535a
    for (i = 0; i < n; ++i) et[i] = iwork[i];
kusano 7d535a
kusano 7d535a
    /* Compute the number of descendants of each node in the etree */
kusano 7d535a
    ifill (relax_end, n, EMPTY);
kusano 7d535a
    ifill (relax_fsupc, n, EMPTY);
kusano 7d535a
    for (j = 0; j < n; j++) descendants[j] = 0;
kusano 7d535a
    for (j = 0; j < n; j++) {
kusano 7d535a
	parent = et[j];
kusano 7d535a
	if ( parent != n )  /* not the dummy root */
kusano 7d535a
	    descendants[parent] += descendants[j] + 1;
kusano 7d535a
    }
kusano 7d535a
kusano 7d535a
    /* Identify the relaxed supernodes by postorder traversal of the etree. */
kusano 7d535a
    for ( f = j = 0; j < n; ) {
kusano 7d535a
	parent = et[j];
kusano 7d535a
	snode_start = j;
kusano 7d535a
	while ( parent != n && descendants[parent] < relax_columns ) {
kusano 7d535a
	    j = parent;
kusano 7d535a
	    parent = et[j];
kusano 7d535a
	}
kusano 7d535a
	/* Found a supernode in postordered etree; j is the last column. */
kusano 7d535a
	++nsuper_et_post;
kusano 7d535a
	k = n;
kusano 7d535a
	for (i = snode_start; i <= j; ++i)
kusano 7d535a
	    k = SUPERLU_MIN(k, inv_post[i]);
kusano 7d535a
	l = inv_post[j];
kusano 7d535a
	if ( (l - k) == (j - snode_start) ) {
kusano 7d535a
	    /* It's also a supernode in the original etree */
kusano 7d535a
	    relax_end[k] = l;		/* Last column is recorded */
kusano 7d535a
	    relax_fsupc[f++] = k;
kusano 7d535a
	    ++nsuper_et;
kusano 7d535a
	} else {
kusano 7d535a
	    for (i = snode_start; i <= j; ++i) {
kusano 7d535a
		l = inv_post[i];
kusano 7d535a
		if ( descendants[i] == 0 ) {
kusano 7d535a
		    relax_end[l] = l;
kusano 7d535a
		    relax_fsupc[f++] = l;
kusano 7d535a
		    ++nsuper_et;
kusano 7d535a
		}
kusano 7d535a
	    }
kusano 7d535a
	}
kusano 7d535a
	j++;
kusano 7d535a
	/* Search for a new leaf */
kusano 7d535a
	while ( descendants[j] != 0 && j < n ) j++;
kusano 7d535a
    }
kusano 7d535a
kusano 7d535a
#if ( PRNTlevel>=1 )
kusano 7d535a
    printf(".. heap_snode_relax:\n"
kusano 7d535a
	   "\tNo of relaxed snodes in postordered etree:\t%d\n"
kusano 7d535a
	   "\tNo of relaxed snodes in original etree:\t%d\n",
kusano 7d535a
	   nsuper_et_post, nsuper_et);
kusano 7d535a
#endif
kusano 7d535a
kusano 7d535a
    /* Recover the original etree */
kusano 7d535a
    for (i = 0; i < n; ++i) et[i] = et_save[i];
kusano 7d535a
kusano 7d535a
    SUPERLU_FREE(post);
kusano 7d535a
    SUPERLU_FREE(iwork);
kusano 7d535a
}