Functions
PORD_INT	greg_pord (PORD_INT, PORD_INT, PORD_INT , PORD_INT , PORD_INT , PORD_INT , PORD_INT *)
graph_t *	newGraph (PORD_INT, PORD_INT)
void	freeGraph (graph_t *)
void	printGraph (graph_t *)
void	randomizeGraph (graph_t *)
graph_t *	setupSubgraph (graph_t , PORD_INT , PORD_INT, PORD_INT *)
graph_t *	setupGraphFromMtx (inputMtx_t *)
graph_t *	setupGridGraph (PORD_INT, PORD_INT, PORD_INT)
PORD_INT	connectedComponents (graph_t *)
graph_t *	compressGraph (graph_t , PORD_INT )
gbisect_t *	newGbisect (graph_t *)
void	freeGbisect (gbisect_t *)
void	printGbisect (gbisect_t *)
void	checkSeparator (gbisect_t *)
void	constructSeparator (gbisect_t , options_t , timings_t *)
PORD_INT	smoothBy2Layers (gbisect_t , PORD_INT , PORD_INT *, PORD_INT, PORD_INT)
void	smoothSeparator (gbisect_t , options_t )
domdec_t *	newDomainDecomposition (PORD_INT, PORD_INT)
void	freeDomainDecomposition (domdec_t *)
void	printDomainDecomposition (domdec_t *)
void	checkDomainDecomposition (domdec_t *)
void	buildInitialDomains (graph_t , PORD_INT , PORD_INT , PORD_INT )
void	mergeMultisecs (graph_t G, PORD_INT , PORD_INT *)
domdec_t *	initialDomainDecomposition (graph_t , PORD_INT , PORD_INT , PORD_INT )
domdec_t *	constructDomainDecomposition (graph_t , PORD_INT )
void	computePriorities (domdec_t , PORD_INT , PORD_INT *, PORD_INT)
void	eliminateMultisecs (domdec_t , PORD_INT , PORD_INT *)
void	findIndMultisecs (domdec_t , PORD_INT , PORD_INT *)
domdec_t *	coarserDomainDecomposition (domdec_t , PORD_INT )
void	shrinkDomainDecomposition (domdec_t *, PORD_INT)
void	checkDDSep (domdec_t *)
PORD_INT	findPseudoPeripheralDomain (domdec_t *, PORD_INT)
void	constructLevelSep (domdec_t *, PORD_INT)
void	initialDDSep (domdec_t *)
void	updateB2W (bucket_t , bucket_t , domdec_t , PORD_INT, PORD_INT , PORD_INT , PORD_INT , PORD_INT *)
void	updateW2B (bucket_t , bucket_t , domdec_t , PORD_INT, PORD_INT , PORD_INT , PORD_INT , PORD_INT *)
void	improveDDSep (domdec_t *)
gbipart_t *	newBipartiteGraph (PORD_INT, PORD_INT, PORD_INT)
void	freeBipartiteGraph (gbipart_t *)
void	printGbipart (gbipart_t *)
gbipart_t *	setupBipartiteGraph (graph_t , PORD_INT , PORD_INT, PORD_INT, PORD_INT *)
void	maximumMatching (gbipart_t , PORD_INT )
void	maximumFlow (gbipart_t , PORD_INT , PORD_INT *)
void	DMviaMatching (gbipart_t , PORD_INT , PORD_INT , PORD_INT )
void	DMviaFlow (gbipart_t , PORD_INT , PORD_INT , PORD_INT , PORD_INT *)
nestdiss_t *	newNDnode (graph_t , PORD_INT , PORD_INT)
void	freeNDnode (nestdiss_t *)
nestdiss_t *	setupNDroot (graph_t , PORD_INT )
void	splitNDnode (nestdiss_t , options_t , timings_t *)
void	buildNDtree (nestdiss_t , options_t , timings_t *)
void	freeNDtree (nestdiss_t *)
multisector_t *	newMultisector (graph_t *)
void	freeMultisector (multisector_t *)
multisector_t *	trivialMultisector (graph_t *)
multisector_t *	constructMultisector (graph_t , options_t , timings_t *)
multisector_t *	extractMS2stage (nestdiss_t *)
multisector_t *	extractMSmultistage (nestdiss_t *)
gelim_t *	newElimGraph (PORD_INT, PORD_INT)
void	freeElimGraph (gelim_t *)
void	printElimGraph (gelim_t *)
gelim_t *	setupElimGraph (graph_t *)
PORD_INT	crunchElimGraph (gelim_t *)
void	buildElement (gelim_t *Gelim, PORD_INT me)
void	updateAdjncy (gelim_t , PORD_INT , PORD_INT, PORD_INT , PORD_INT )
void	findIndNodes (gelim_t , PORD_INT , PORD_INT, PORD_INT , PORD_INT , PORD_INT , PORD_INT )
void	updateDegree (gelim_t , PORD_INT , PORD_INT, PORD_INT *)
void	updateScore (gelim_t , PORD_INT , PORD_INT, PORD_INT, PORD_INT *)
elimtree_t *	extractElimTree (gelim_t *)
bucket_t *	newBucket (PORD_INT, PORD_INT, PORD_INT)
void	freeBucket (bucket_t *)
bucket_t *	setupBucket (PORD_INT, PORD_INT, PORD_INT)
PORD_INT	minBucket (bucket_t *)
void	insertBucket (bucket_t *, PORD_INT, PORD_INT)
void	removeBucket (bucket_t *, PORD_INT)
minprior_t *	newMinPriority (PORD_INT nvtx, PORD_INT nstages)
void	freeMinPriority (minprior_t *)
minprior_t *	setupMinPriority (multisector_t *)
elimtree_t *	orderMinPriority (minprior_t , options_t , timings_t *)
void	eliminateStage (minprior_t , PORD_INT, PORD_INT, timings_t )
PORD_INT	eliminateStep (minprior_t *, PORD_INT, PORD_INT)
elimtree_t *	newElimTree (PORD_INT, PORD_INT)
void	freeElimTree (elimtree_t *)
void	printElimTree (elimtree_t *)
PORD_INT	firstPostorder (elimtree_t *)
PORD_INT	firstPostorder2 (elimtree_t *, PORD_INT)
PORD_INT	nextPostorder (elimtree_t *, PORD_INT)
PORD_INT	firstPreorder (elimtree_t *)
PORD_INT	nextPreorder (elimtree_t *, PORD_INT)
elimtree_t *	setupElimTree (graph_t , PORD_INT , PORD_INT *)
void	initFchSilbRoot (elimtree_t *)
void	permFromElimTree (elimtree_t , PORD_INT )
elimtree_t *	expandElimTree (elimtree_t , PORD_INT , PORD_INT)
elimtree_t *	permuteElimTree (elimtree_t , PORD_INT )
elimtree_t *	fundamentalFronts (elimtree_t *)
elimtree_t *	mergeFronts (elimtree_t *, PORD_INT)
elimtree_t *	compressElimTree (elimtree_t , PORD_INT , PORD_INT)
PORD_INT	justifyFronts (elimtree_t *)
PORD_INT	nWorkspace (elimtree_t *)
PORD_INT	nFactorIndices (elimtree_t *)
PORD_INT	nFactorEntries (elimtree_t *)
FLOAT	nFactorOps (elimtree_t *)
void	subtreeFactorOps (elimtree_t , FLOAT )
FLOAT	nTriangularOps (elimtree_t *)
inputMtx_t *	newInputMtx (PORD_INT, PORD_INT)
void	freeInputMtx (inputMtx_t *)
void	printInputMtx (inputMtx_t *)
denseMtx_t *	newDenseMtx (workspace_t *, PORD_INT)
void	freeDenseMtx (denseMtx_t *)
void	printDenseMtx (denseMtx_t *)
void	checkDenseMtx (denseMtx_t *)
workspace_t *	initWorkspaceForDenseMtx (PORD_INT, PORD_INT)
FLOAT *	getWorkspaceForDenseMtx (workspace_t *, PORD_INT)
void	freeWorkspaceForDenseMtx (workspace_t *)
inputMtx_t *	setupInputMtxFromGraph (graph_t *)
inputMtx_t *	setupLaplaceMtx (PORD_INT, PORD_INT, PORD_INT)
inputMtx_t *	permuteInputMtx (inputMtx_t , PORD_INT )
css_t *	newCSS (PORD_INT, PORD_INT, PORD_INT)
void	freeCSS (css_t *)
css_t *	setupCSSFromGraph (graph_t , PORD_INT , PORD_INT *)
css_t *	setupCSSFromFrontSubscripts (frontsub_t *)
frontsub_t *	newFrontSubscripts (elimtree_t *)
void	freeFrontSubscripts (frontsub_t *)
void	printFrontSubscripts (frontsub_t *)
frontsub_t *	setupFrontSubscripts (elimtree_t , inputMtx_t )
factorMtx_t *	newFactorMtx (PORD_INT)
void	freeFactorMtx (factorMtx_t *)
void	printFactorMtx (factorMtx_t *)
void	initFactorMtx (factorMtx_t L, inputMtx_t )
void	initFactorMtxNEW (factorMtx_t L, inputMtx_t )
void	numfac (factorMtx_t L, timings_t cpus)
denseMtx_t *	setupFrontalMtx (workspace_t , factorMtx_t , PORD_INT)
void	initLocalIndices (denseMtx_t , PORD_INT , PORD_INT *)
denseMtx_t *	extendedAdd (denseMtx_t , denseMtx_t , PORD_INT , PORD_INT )
denseMtx_t *	setupUpdateMtxFromFrontalMtx (denseMtx_t , factorMtx_t )
denseMtx_t *	factorize1x1Kernel (denseMtx_t *, PORD_INT)
denseMtx_t *	factorize2x2Kernel (denseMtx_t *, PORD_INT)
denseMtx_t *	factorize3x3Kernel (denseMtx_t *, PORD_INT)
void	forwardSubst1x1 (factorMtx_t , FLOAT )
void	backwardSubst1x1 (factorMtx_t , FLOAT )
void	forwardSubst1x1NEW (factorMtx_t , FLOAT )
void	backwardSubst1x1NEW (factorMtx_t , FLOAT )
mapping_t *	newMapping (elimtree_t *, PORD_INT)
void	freeMapping (mapping_t *)
void	printMapping (mapping_t *)
void	listing (mapping_t , PORD_INT, PORD_INT, PORD_INT, FLOAT , FLOAT *)
mapping_t *	setupMapping (elimtree_t *, PORD_INT, PORD_INT)
void	split (mapping_t , PORD_INT, PORD_INT, PORD_INT, PORD_INT , PORD_INT , FLOAT , PORD_INT)
elimtree_t *	SPACE_ordering (graph_t , options_t , timings_t *)
elimtree_t *	SPACE_transformElimTree (elimtree_t *, PORD_INT)
factorMtx_t *	SPACE_symbFac (elimtree_t , inputMtx_t )
void	SPACE_numFac (factorMtx_t , timings_t )
void	SPACE_solveTriangular (factorMtx_t L, FLOAT rhs, FLOAT *xvec)
void	SPACE_solve (inputMtx_t , FLOAT , FLOAT , options_t , timings_t *)
void	SPACE_solveWithPerm (inputMtx_t , PORD_INT , FLOAT , FLOAT , options_t , timings_t )
mapping_t *	SPACE_mapping (graph_t , PORD_INT , options_t , timings_t )
void	insertUpInts (PORD_INT, PORD_INT *)
void	insertUpIntsWithStaticIntKeys (PORD_INT, PORD_INT , PORD_INT )
void	insertDownIntsWithStaticFloatKeys (PORD_INT, PORD_INT , FLOAT )
void	insertUpFloatsWithIntKeys (PORD_INT, FLOAT , PORD_INT )
void	qsortUpInts (PORD_INT, PORD_INT , PORD_INT )
void	qsortUpFloatsWithIntKeys (PORD_INT, FLOAT , PORD_INT , PORD_INT *)
void	distributionCounting (PORD_INT, PORD_INT , PORD_INT )
graph_t *	readChacoGraph (char *)
inputMtx_t *	readHarwellBoeingMtx (char *)
topology_t *	newTopology (PORD_INT)
void	freeTopology (topology_t *)
void	printTopology (topology_t *)
topology_t *	setupTopology (void)
void	recMapCube (topology_t *, PORD_INT, PORD_INT, PORD_INT, PORD_INT, PORD_INT, PORD_INT)
void	sendCube (topology_t , void , size_t, PORD_INT)
size_t	recvCube (topology_t , void , size_t, PORD_INT)
PORD_INT	myrank (void)
mask_t *	newMask (PORD_INT)
void	freeMask (mask_t *)
mask_t *	setupMask (PORD_INT, PORD_INT, PORD_INT)
void	broadcastInputMtx (topology_t , inputMtx_t *)
void	broadcastElimTree (topology_t , elimtree_t *)
void	broadcastArray (topology_t , char , size_t)
buffer_t *	newBuffer (size_t)
void	freeBuffer (buffer_t *)
buffer_t *	exchangeBuffer (topology_t , buffer_t , PORD_INT)
buffer_t *	setupSymbFacBuffer (frontsub_t , PORD_INT )
void	readoutSymbFacBuffer (buffer_t , frontsub_t , PORD_INT *)
buffer_t *	setupNumFacBuffer (workspace_t , mask_t , PORD_INT)
void	readoutNumFacBuffer (workspace_t , buffer_t , denseMtx_t **)
buffer_t *	setupTriangularBuffer (frontsub_t , PORD_INT , FLOAT *)
void	readoutTriangularBuffer (buffer_t , frontsub_t , PORD_INT , FLOAT )
frontsub_t *	newFrontSubscriptsPAR (mask_t , mapping_t , elimtree_t *)
frontsub_t *	setupFrontSubscriptsPAR (topology_t , mask_t , mapping_t , elimtree_t , inputMtx_t *)
css_t *	setupCSSFromFrontSubscriptsPAR (mask_t , mapping_t , frontsub_t *)
void	initFactorMtxPAR (mask_t , mapping_t , factorMtx_t , inputMtx_t )
void	numfacPAR (topology_t , mask_t , mapping_t , factorMtx_t , PORD_INT msglvl, timings_t *)
denseMtx_t *	setupFrontalMtxPAR (mask_t , PORD_INT, workspace_t , factorMtx_t *, PORD_INT)
void	initLocalIndicesPAR (denseMtx_t , PORD_INT , PORD_INT *)
denseMtx_t *	extendedAddPAR (denseMtx_t , denseMtx_t , PORD_INT , PORD_INT )
denseMtx_t *	setupUpdateMtxFromFrontalMtxPAR (denseMtx_t , factorMtx_t )
denseMtx_t *	setupUpdateMtxFromBuffer (workspace_t , FLOAT )
void	splitDenseMtxColumnWise (denseMtx_t , mask_t , buffer_t *, PORD_INT)
void	splitDenseMtxRowWise (denseMtx_t , mask_t , buffer_t *, PORD_INT)
denseMtx_t *	factorize1x1KernelPAR (topology_t , mask_t , PORD_INT, denseMtx_t , frontsub_t , timings_t *)
denseMtx_t *	factorize2x2KernelPAR (topology_t , mask_t , PORD_INT, denseMtx_t , frontsub_t , timings_t *)
denseMtx_t *	factorize3x3KernelPAR (topology_t , mask_t , PORD_INT, denseMtx_t , frontsub_t , timings_t *)
void	forwardSubst1x1PAR (topology_t , mask_t , mapping_t , factorMtx_t , FLOAT , FLOAT )
void	backwardSubst1x1PAR (topology_t , mask_t , mapping_t , factorMtx_t , FLOAT *)
void	forwardSubst1x1KernelPAR (topology_t , mask_t , PORD_INT, PORD_INT, factorMtx_t , FLOAT , FLOAT *)
void	backwardSubst1x1KernelPAR (topology_t , mask_t , PORD_INT, PORD_INT, factorMtx_t , FLOAT )
void	accumulateVector (topology_t , mask_t , mapping_t , factorMtx_t , FLOAT *)
topology_t *	SPACE_setupTopology (void)
mask_t *	SPACE_setupMask (topology_t *, PORD_INT)
void	SPACE_cleanup (topology_t , mask_t )
factorMtx_t *	SPACE_symbFacPAR (topology_t , mask_t , mapping_t , elimtree_t , inputMtx_t *)
void	SPACE_numFacPAR (topology_t , mask_t , mapping_t , factorMtx_t , PORD_INT msglvl, timings_t *)
void	SPACE_solveTriangularPAR (topology_t , mask_t , mapping_t , factorMtx_t , FLOAT , FLOAT )
void	SPACE_solveWithPermPAR (topology_t top, mask_t mask, inputMtx_t A, PORD_INT perm, FLOAT rhs, FLOAT xvec, options_t options, timings_t cpus)

Function Documentation

◆ accumulateVector()

void accumulateVector	(	topology_t *	,
		mask_t *	,
		mapping_t *	,
		factorMtx_t *	,
		FLOAT *	)

◆ backwardSubst1x1()

void backwardSubst1x1	(	factorMtx_t *	,
		FLOAT *	)

◆ backwardSubst1x1KernelPAR()

void backwardSubst1x1KernelPAR	(	topology_t *	,
		mask_t *	,
		PORD_INT	,
		PORD_INT	,
		factorMtx_t *	,
		FLOAT *	)

◆ backwardSubst1x1NEW()

void backwardSubst1x1NEW	(	factorMtx_t *	,
		FLOAT *	)

◆ backwardSubst1x1PAR()

void backwardSubst1x1PAR	(	topology_t *	,
		mask_t *	,
		mapping_t *	,
		factorMtx_t *	,
		FLOAT *	)

◆ broadcastArray()

void broadcastArray	(	topology_t *	,
		char *	,
		size_t	)

◆ broadcastElimTree()

void broadcastElimTree	(	topology_t *	,
		elimtree_t **	)

◆ broadcastInputMtx()

void broadcastInputMtx	(	topology_t *	,
		inputMtx_t **	)

◆ buildElement()

void buildElement	(	gelim_t *	Gelim,
		PORD_INT	me )

Definition at line 357 of file gelim.c.

{ graph_t *G;
  PORD_INT     *xadj, *adjncy, *vwght, *len, *elen, *parent, *degree, *score;
  PORD_INT     degme, elenme, vlenme, mesrcptr, medeststart, medeststart2;
  PORD_INT     medestptr, ln, p, i, j, v, e;
 
  G = Gelim->G;
  xadj = G->xadj;
  adjncy = G->adjncy;
  vwght = G->vwght;
  len = Gelim->len;
  elen = Gelim->elen;
  parent = Gelim->parent;
  degree = Gelim->degree;
  score = Gelim->score;
 
  /* ---------------------------------
     construct boundary of element Lme
     --------------------------------- */
  degme = 0;
  G->totvwght -= vwght[me];  /* me eliminated => reduce weight of Gelim */
  vwght[me] = -vwght[me];
  score[me] = -3;            /* variable me becomes an element */
 
  elenme = elen[me];
  vlenme = len[me] - elenme;
  mesrcptr = xadj[me];
 
  /* -----------------------------------------------------------
     if me is a leaf => its boundary can be constructed in-place
     ----------------------------------------------------------- */
  if (elenme == 0)
   { medeststart = xadj[me];         /* Lme overwrites old variable */
     medestptr = medeststart;        /* boundary of Lme starts here */
     for (i = 0; i < vlenme; i++)
      { v = adjncy[mesrcptr++];
        if (vwght[v] > 0)            /* v not yet placed in boundary */
         { degme += vwght[v];        /* increase size of Lme */
           vwght[v] = -vwght[v];     /* flag v as being in Lme */
           adjncy[medestptr++] = v;
         }
      }
   }
 
  /* -------------------------------------------------------------------
     me is not a leaf => its boundary must be constructed in empty space
     ------------------------------------------------------------------- */
  else
   { medeststart = G->nedges;        /* Lme appended to graph */
     medestptr = medeststart;        /* boundary of Lme starts here */
     for (i = 0; i <= elenme; i++)
      { if (i < elenme)              /* working on elements */
         { len[me]--;
           e = adjncy[mesrcptr++];   /* merge boundary of element e with Lme */
           p = xadj[e];              /* adjacency list of e starts here */
           ln = len[e];
         }
        else
         { e = me;                   /* merge uncovered variables with Lme */
           p = mesrcptr;             /* variables start here */
           ln = vlenme;
         }
        for (j = 0; j < ln; j++)
         { len[e]--;                 /* pick next variable, decrease length */
           v = adjncy[p++];
           if (vwght[v] > 0)
            { degme += vwght[v];     /* increase size of Lme */
              vwght[v] = -vwght[v];  /* flag v as being in Lme */
 
              /* ------------------------------------------
                 add v to Lme, compress adjncy if necessary
                 ------------------------------------------ */
              if (medestptr == Gelim->maxedges)
               { if (len[me] == 0) xadj[me] = -1;
                 else xadj[me] = mesrcptr;
                 if (len[e] == 0) xadj[e] = -1;
                 else xadj[e] = p;
 
                 /* crunch adjacency list -- !!!we need more memory!!! */
                 if (!crunchElimGraph(Gelim))
                  { fprintf(stderr, "\nError in function buildElement\n"
                         "  unable to construct element (not enough memory)\n");
                    quit();
                  }
 
                 /* crunch partially constructed element me */
                 medeststart2 = G->nedges;
                 for (p = medeststart; p < medestptr; p++)
                   adjncy[G->nedges++] = adjncy[p];
                 medeststart = medeststart2;
                 medestptr = G->nedges;
 
                 mesrcptr = xadj[me];
                 p = xadj[e];
               }
              adjncy[medestptr++] = v;
            }
         }
 
        /* ----------------------
           mark absorbed elements
           ---------------------- */
        if (e != me)
         { xadj[e] = -1;
           parent[e] = me;
           score[e] = -4;
         }
      }
 
     G->nedges = medestptr;          /* new element Lme ends here */
   }
 
  /* -----------------------------------
     element me successfully constructed
     ----------------------------------- */
  degree[me] = degme;
  xadj[me] = medeststart;
  vwght[me] = -vwght[me];
  elen[me] = 0;
  len[me] = medestptr - medeststart;
  if (len[me] == 0)
    xadj[me] = -1;
 
  /* ---------------------------
     unmark all variables in Lme
     --------------------------- */
  mesrcptr = xadj[me];
  vlenme = len[me];
  for (i = 0; i < vlenme; i++)
   { v = adjncy[mesrcptr++];
     vwght[v] = -vwght[v];
   }
}

◆ buildInitialDomains()

void buildInitialDomains	(	graph_t *	G,
		PORD_INT *	vtxlist,
		PORD_INT *	vtype,
		PORD_INT *	rep )

Definition at line 224 of file ddcreate.c.

{ PORD_INT *xadj, *adjncy;
  PORD_INT nvtx, u, v, w, i, j, jstart, jstop;
 
  xadj = G->xadj;
  adjncy = G->adjncy;
  nvtx = G->nvtx;
 
  /* --------------------------------------------------------------------
     determine initial domains according to the order of nodes in vtxlist
     -------------------------------------------------------------------- */
  for (i = 0; i < nvtx; i++)
   { u = vtxlist[i];
     if (vtype[u] == 0)
      { vtype[u] = 1;
        jstart = xadj[u];
        jstop = xadj[u+1];
        for (j = jstart; j < jstop; j++)
         { v = adjncy[j];
           vtype[v] = 2;
         }
      }
   }
 
  /* ------------------------------------------------------------
     eliminate all multisecs that are adjacent to only one domain
     ------------------------------------------------------------ */
  for (i = 0; i < nvtx; i++)
   { u = vtxlist[i];
     if (vtype[u] == 2)
      { v = -1;
        jstart = xadj[u];
        jstop = xadj[u+1];
        for (j = jstart; j < jstop; j++)
         { w = adjncy[j];
           if (vtype[w] == 1)
            { if (v == -1)
                v = rep[w];          /* u adjacent to domain v = rep[w] */
              else if (v != rep[w])
               { v = -1;             /* u adjacent to another domain */
                 break;
               }
            }
         }
        if (v != -1)                 /* u absorbed by domain v */
         { vtype[u] = 1;
           rep[u] = v;
         }
      }
   }
}

◆ buildNDtree()

void buildNDtree	(	nestdiss_t *	ndroot,
		options_t *	options,
		timings_t *	cpus )

Definition at line 213 of file nestdiss.c.

{ nestdiss_t *nd;
  nestdiss_t *queue[2*MAX_SEPS+1];
  PORD_INT        maxseps, seps, domainsize, qhead, qtail;
 
  maxseps = MAX_SEPS;
  domainsize = options[OPTION_DOMAIN_SIZE];
  if (domainsize == 1) maxseps = DEFAULT_SEPS;   /* secret switch */
 
  /* --------------------------------------------------
     build the nested dissection tree under root ndroot
     -------------------------------------------------- */
  queue[0] = ndroot;
  qhead = 0; qtail = 1; seps = 0;
  while ((qhead != qtail) && (seps < maxseps))
   { seps++;
     nd = queue[qhead++];
 
     splitNDnode(nd, options, cpus);
     if ((nd->childB == NULL) || (nd->childW == NULL))
      { fprintf(stderr, "\nError in function buildNDtree\n"
             "  recursive nested dissection process failed\n");
        quit();
      }
 
     if (options[OPTION_MSGLVL] > 1)
       printf("%4d. S %6d, B %6d, W %6d [bal %4.2f, rel %6.4f, cost %7.2f]\n",
             seps, nd->cwght[GRAY], nd->cwght[BLACK], nd->cwght[WHITE],
             (FLOAT)min(nd->cwght[BLACK], nd->cwght[WHITE])
                     / max(nd->cwght[BLACK], nd->cwght[WHITE]),
             (FLOAT)nd->cwght[GRAY]
                     / (nd->cwght[GRAY] + nd->cwght[BLACK] + nd->cwght[WHITE]),
             F(nd->cwght[GRAY], nd->cwght[BLACK], nd->cwght[WHITE]));
 
     if ((nd->childB->nvint > MIN_NODES)
        && ((nd->cwght[BLACK] > domainsize) || (qtail < DEFAULT_SEPS)))
       queue[qtail++] = nd->childB;
     if ((nd->childW->nvint > MIN_NODES)
        && ((nd->cwght[WHITE] > domainsize) || (qtail < DEFAULT_SEPS)))
       queue[qtail++] = nd->childW;
   }
}

◆ checkDDSep()

void checkDDSep ( domdec_t * dd )

Definition at line 68 of file ddbisect.c.

{ PORD_INT *xadj, *adjncy, *vwght, *vtype, *color, *cwght;
  PORD_INT nvtx, err, u, v, i, istart, istop, nBdom, nWdom;
  PORD_INT checkS, checkB, checkW;
 
  nvtx = dd->G->nvtx;
  xadj = dd->G->xadj;
  adjncy = dd->G->adjncy;
  vwght = dd->G->vwght;
  vtype = dd->vtype;
  color = dd->color;
  cwght = dd->cwght;
 
  err = FALSE;
  printf("checking separator of domain decomposition (S %d, B %d, W %d)\n",
         cwght[GRAY], cwght[BLACK], cwght[WHITE]);
 
  checkS = checkB = checkW = 0;
  for (u = 0; u < nvtx; u++)
    /* check neighborhood of multisector nodes */
    if (vtype[u] == 2)
     { nBdom = nWdom = 0;
       istart = xadj[u];
       istop = xadj[u+1];
       for (i = istart; i < istop; i++)
        { v = adjncy[i];
          if (color[v] == BLACK) nBdom++;
          if (color[v] == WHITE) nWdom++;
        }
       switch(color[u])
        { case GRAY:
            checkS += vwght[u];
            if ((nBdom == 0) || (nWdom == 0))
              printf("WARNING: multisec %d belongs to S, but nBdom = %d and "
                     "nWdom = %d\n", u, nBdom, nWdom);
            break;
          case BLACK: 
            checkB += vwght[u];
            if (nWdom > 0)
             { printf("ERROR: black multisec %d adjacent to white domain\n", u);
               err = TRUE;
             }
            break;
          case WHITE:
            checkW += vwght[u];
            if (nBdom > 0)
             { printf("ERROR: white multisec %d adjacent to black domain\n", u);
               err = TRUE;
             }
            break;
          default:
            printf("ERROR: multisec %d has unrecognized color %d\n", u,
                   color[u]);
            err = TRUE;
        }
     }
 
    /* sum up size of white/black domains */
    else /* if (vtype[u] == 1) */
      switch(color[u])
       { case BLACK:
           checkB += vwght[u]; break;
         case WHITE:
           checkW += vwght[u]; break;
         default:
           printf("ERROR: domain %d has unrecognized color %d\n", u, color[u]);
           err = TRUE;
       }
 
  /* check cwght[GRAY], cwght[BLACK], cwght[WHITE] */
  if ((checkS != cwght[GRAY]) || (checkB != cwght[BLACK])
     || (checkW != cwght[WHITE]))
   { printf("ERROR in partitioning: checkS %d (S %d), checkB %d (B %d), "
            "checkW %d (W %d)\n", checkS, cwght[GRAY], checkB, cwght[BLACK],
             checkW, cwght[WHITE]);
     err = TRUE;
   }
  if (err) quit();
}

◆ checkDenseMtx()

void checkDenseMtx ( denseMtx_t * )

◆ checkDomainDecomposition()

void checkDomainDecomposition ( domdec_t * dd )

Definition at line 163 of file ddcreate.c.

{ PORD_INT *xadj, *adjncy, *vwght, *vtype;
  PORD_INT err, nvtx, ndom, domwght, dom, multi, u, v, i, istart, istop;
 
  nvtx = dd->G->nvtx;
  xadj = dd->G->xadj;
  adjncy = dd->G->adjncy;
  vwght = dd->G->vwght;
  vtype = dd->vtype;
 
  err = FALSE;
  printf("checking domain decomposition (#nodes %d, #edges %d)\n",
         dd->G->nvtx, dd->G->nedges >> 1);
 
  ndom = domwght = 0;
  for (u = 0; u < nvtx; u++)
   { /* check node type */
     if ((vtype[u] != 1) && (vtype[u] != 2))
      { printf("ERROR: node %d is neither DOMAIN nor MULTISEC\n", u);
        err = TRUE;
      }
     /* count domains and sum up their weight */
     if (vtype[u] == 1)
      { ndom++;
        domwght += vwght[u];
      } 
     /* check number of neighboring domains and multisecs */
     dom = multi = 0;
     istart = xadj[u];
     istop = xadj[u+1];
     for (i = istart; i < istop; i++)
      { v = adjncy[i];
        if (vtype[v] == 1) dom++;
        if (vtype[v] == 2) multi++;
      }
     if ((vtype[u] == 1) && (dom > 0))
      { printf("ERROR: domain %d is adjacent to other domain\n", u);
        err = TRUE;
      }
     if ((vtype[u] == 2) && (dom < 2))
      { printf("ERROR: less than 2 domains adjacent to multisec node %d\n", u);
        err = TRUE;
      }
     if ((vtype[u] == 2) && (multi > 0))
      { printf("ERROR: multisec %d is adjacent to other multisec nodes\n", u);
        err = TRUE;
      }
   }
  /* check number and weight of domains */
  if ((ndom != dd->ndom) || (domwght != dd->domwght))
   { printf("ERROR: number/size (%d/%d) of domains does not match with those in"
            " domain decomp. (%d/%d)\n", ndom, domwght, dd->ndom, dd->domwght);
     err = TRUE;
   }
  if (err) quit();
}

◆ checkSeparator()

void checkSeparator ( gbisect_t * Gbisect )

Definition at line 117 of file gbisect.c.

{ PORD_INT *xadj, *adjncy, *vwght, *color, *cwght;
  PORD_INT nvtx, err, checkS, checkB, checkW, a, b, u, v, i, istart, istop;
 
  nvtx = Gbisect->G->nvtx;
  xadj = Gbisect->G->xadj;
  adjncy = Gbisect->G->adjncy;
  vwght = Gbisect->G->vwght;
  color = Gbisect->color;
  cwght = Gbisect->cwght;
 
  err = FALSE;
  printf("checking separator of induced subgraph (S %d, B %d, W %d)\n",
         cwght[GRAY], cwght[BLACK], cwght[WHITE]);
 
  checkS = checkB = checkW = 0;
  for (u = 0; u < nvtx; u++)
   { istart = xadj[u];
     istop = xadj[u+1];
     switch(color[u])
      { case GRAY:                /* is it a minimal separator? */
          checkS += vwght[u];
          a = b = FALSE;
          for (i = istart; i < istop; i++)
           { v = adjncy[i];
             if (color[v] == WHITE) a = TRUE;
             if (color[v] == BLACK) b = TRUE;
           }
          if (!((a) && (b)))
            printf("WARNING: not a minimal separator (node %d)\n", u);
          break;
        case BLACK:               /* is it realy a separator? */
          checkB += vwght[u];
          for (i = istart; i < istop; i++)
           { v = adjncy[i];
             if (color[v] == WHITE)
              { printf("ERROR: white node %d adjacent to black node %d\n", u,v);
                err = TRUE;
              }
           }
          break;
        case WHITE:
          checkW += vwght[u];
          break;
        default:
          printf("ERROR: node %d has unrecognized color %d\n", u, color[u]);
          err = TRUE;
      }
   }
 
  /* check cwght[GRAY], cwght[BLACK], cwght[WHITE] */
  if ((checkS != cwght[GRAY]) || (checkB != cwght[BLACK])
     || (checkW != cwght[WHITE]))
   { printf("ERROR in partitioning: checkS %d (S %d), checkB %d (B %d), "
            "checkW %d (W %d)\n", checkS, cwght[GRAY], checkB, cwght[BLACK],
             checkW, cwght[WHITE]);
     err = TRUE;
   }
  if (err) quit();
}

◆ coarserDomainDecomposition()

domdec_t * coarserDomainDecomposition	(	domdec_t *	dd1,
		PORD_INT *	rep )

Definition at line 778 of file ddcreate.c.

{ domdec_t *dd2;
  PORD_INT      *xadjdd1, *adjncydd1, *vwghtdd1, *vtypedd1, *mapdd1;
  PORD_INT      *xadjdd2, *adjncydd2, *vwghtdd2, *vtypedd2;
  PORD_INT      *tmp, *bin, nvtxdd1, nedgesdd1, nvtxdd2, nedgesdd2;
  PORD_INT      ndom, domwght, flag, u, v, w, i, istart, istop;
 
  nvtxdd1 = dd1->G->nvtx;
  nedgesdd1 = dd1->G->nedges;
  xadjdd1 = dd1->G->xadj;
  adjncydd1 = dd1->G->adjncy;
  vwghtdd1 = dd1->G->vwght;
  vtypedd1 = dd1->vtype;
  mapdd1 = dd1->map;
 
  /* ------------------------
     allocate working storage
     ------------------------ */
  mymalloc(tmp, nvtxdd1, PORD_INT);
  mymalloc(bin, nvtxdd1, PORD_INT);
  for (u = 0; u < nvtxdd1; u++)
   { tmp[u] = -1;
     bin[u] = -1;
   }
 
  /* ------------------------------------------------------------
     allocate memory using the upper bounds nvtxdd1 and nedgesdd1
     ------------------------------------------------------------ */
  dd2 = newDomainDecomposition(nvtxdd1, nedgesdd1);
  xadjdd2 = dd2->G->xadj;
  adjncydd2 = dd2->G->adjncy;
  vwghtdd2 = dd2->G->vwght;
  vtypedd2 = dd2->vtype;
 
  /* -------------------------------------------------------
     put all nodes u belonging to representative v in bin[v]
     ------------------------------------------------------- */
  for (u = 0; u < nvtxdd1; u++)
   { v = rep[u];
     if (u != v)
      { bin[u] = bin[v];
        bin[v] = u;
      }
   }
 
  /* ----------------------------------------------
     and now build the coarser domain decomposition
     ---------------------------------------------- */
  flag = 1;
  nvtxdd2 = nedgesdd2 = 0;
  ndom = domwght = 0;
  for (u = 0; u < nvtxdd1; u++)
    if (rep[u] == u)
     { xadjdd2[nvtxdd2] = nedgesdd2;
       vwghtdd2[nvtxdd2] = 0;
       vtypedd2[nvtxdd2] = vtypedd1[u];
       if (vtypedd2[nvtxdd2] == 3)
         vtypedd2[nvtxdd2] = 1;
       tmp[u] = flag;
 
       /* find all cluster that are adjacent to u in dom. dec. */
       v = u;
       do
        { mapdd1[v] = nvtxdd2;
          vwghtdd2[nvtxdd2] += vwghtdd1[v];
          if ((vtypedd1[v] == 1) || (vtypedd1[v] == 2))
           { istart = xadjdd1[v];
             istop = xadjdd1[v+1];
             for (i = istart; i < istop; i++)
              { w = adjncydd1[i];
                if (tmp[rep[w]] != flag)
                 { tmp[rep[w]] = flag;
                   adjncydd2[nedgesdd2++] = rep[w];
                 }
              }
           }
          v = bin[v];
        } while (v != -1);
 
       if (vtypedd2[nvtxdd2] == 1)
        { ndom++;
          domwght += vwghtdd2[nvtxdd2];
        }
       nvtxdd2++;
       flag++;
     }
 
  /* --------------------------------------------
     finalize the new domain decomposition object
     -------------------------------------------- */
  xadjdd2[nvtxdd2] = nedgesdd2;
  dd2->G->nvtx = nvtxdd2;
  dd2->G->nedges = nedgesdd2;
  dd2->G->type = WEIGHTED;
  dd2->G->totvwght = dd1->G->totvwght;
  for (i = 0; i < nedgesdd2; i++)
    adjncydd2[i] = mapdd1[adjncydd2[i]];
  for (u = 0; u < nvtxdd2; u++)
    dd2->color[u] = dd2->map[u] = -1;
  dd2->ndom = ndom;
  dd2->domwght = domwght;
 
  /* --------------------------
     set back node types in dd1
     -------------------------- */
  for (u = 0; u < nvtxdd1; u++)
    if ((vtypedd1[u] == 3) || (vtypedd1[u] == 4))
      vtypedd1[u] = 2;
 
  /* -------------------------------
     free working storage and return
     ------------------------------- */
  free(tmp); free(bin);
  return(dd2);
}

◆ compressElimTree()

elimtree_t * compressElimTree	(	elimtree_t *	T,
		PORD_INT *	frontmap,
		PORD_INT	cnfronts )

Definition at line 689 of file tree.c.

{ elimtree_t *T2;
  PORD_INT        *ncolfactor, *ncolupdate, *parent, *vtx2front;
  PORD_INT        nvtx, nfronts, u, K, pK, newfront, pnewfront;
 
  nvtx = T->nvtx;
  nfronts = T->nfronts;
  ncolfactor = T->ncolfactor;
  ncolupdate = T->ncolupdate;
  parent = T->parent;
  vtx2front = T->vtx2front;
 
  /* --------------------------------------------
     allocate memory for the new elimtree T2
     and init. ncolfactor, ncolupdate, and parent
     -------------------------------------------- */
  T2 = newElimTree(nvtx, cnfronts);
  for (K = 0; K < cnfronts; K++)
   { T2->ncolfactor[K] = T2->ncolupdate[K] = 0;
     T2->parent[K] = -1;
   }
 
  /* --------------------------------------------------------------
     set the new vectors T2->ncolfactor, T2->ncolupdate, T2->parent
     -------------------------------------------------------------- */
  for (K = 0; K < nfronts; K++)
   { newfront = frontmap[K];
     T2->ncolfactor[newfront] += ncolfactor[K];
     if (((pK = parent[K]) != -1)
        && ((pnewfront = frontmap[pK]) != newfront))
      { T2->parent[newfront] = pnewfront;
        T2->ncolupdate[newfront] = ncolupdate[K];
      }
   }
 
  /* ---------------------------------------------------
     set the new vectors T2->firstchild and T2->silbings
     --------------------------------------------------- */
  initFchSilbRoot(T2);
 
  /* ------------------------------------
     set the the new vector T2->vtx2front
     ------------------------------------ */
  for (u = 0; u < nvtx; u++)
    T2->vtx2front[u] = frontmap[vtx2front[u]];
  return(T2);
}

◆ compressGraph()

graph_t * compressGraph	(	graph_t *	G,
		PORD_INT *	vtxmap )

Definition at line 470 of file graph.c.

{ graph_t *Gc;
  PORD_INT     *xadj, *adjncy, *vwght, *xadjGc, *adjncyGc, *vwghtGc, *perm;
  PORD_INT     nvtx, nvtxGc, nedgesGc, u, v, i, istart, istop;
 
  nvtx = G->nvtx;
  xadj = G->xadj;
  adjncy = G->adjncy;
  vwght = G->vwght;
 
  /* --------------------------------------------------------------
     compressed graph small enough? if so, allocate working storage
     -------------------------------------------------------------- */
  /* avoid print statement
   * printf("indNodes(G, vtxmap) = %d",indNodes(G, vtxmap)); */
  if ((nvtxGc = indNodes(G, vtxmap)) > COMPRESS_FRACTION * nvtx)
    return(NULL);
  mymalloc(perm, nvtx, PORD_INT);
 
  /* -----------------------------------
     count edges of the compressed graph
     ----------------------------------- */
  nedgesGc = 0;
  for (u = 0; u < nvtx; u++)
   if (vtxmap[u] == u)
    { istart = xadj[u];
      istop = xadj[u+1];
      for (i = istart; i < istop; i++)
       { v = adjncy[i];
         if (vtxmap[v] == v) nedgesGc++;
       }
    }
 
  /* ---------------------------------------------------------
     allocate memory for the compressed graph and construct it
     --------------------------------------------------------- */
  Gc = newGraph(nvtxGc, nedgesGc);
  xadjGc = Gc->xadj;
  adjncyGc = Gc->adjncy;
  vwghtGc = Gc->vwght;
 
  nvtxGc = nedgesGc = 0;
  for (u = 0; u < nvtx; u++)
   if (vtxmap[u] == u)
    { istart = xadj[u];
      istop = xadj[u+1];
      xadjGc[nvtxGc] = nedgesGc;
      vwghtGc[nvtxGc] = 0;
      perm[u] = nvtxGc++;
      for (i = istart; i < istop; i++)
       { v = adjncy[i];
         if (vtxmap[v] == v) adjncyGc[nedgesGc++] = v;
       }
    }
  xadjGc[nvtxGc] = nedgesGc;
 
  for (i = 0; i < nedgesGc; i++)
    adjncyGc[i] = perm[adjncyGc[i]];
  for (u = 0; u < nvtx; u++)
   { vtxmap[u] = perm[vtxmap[u]];
     vwghtGc[vtxmap[u]] += vwght[u];
   }
  Gc->type = WEIGHTED;
  Gc->totvwght = G->totvwght;
 
  /* ----------------------
     free memory and return
     ---------------------- */
  free(perm);
  return(Gc);
}

◆ computePriorities()

void computePriorities	(	domdec_t *	dd,
		PORD_INT *	msvtxlist,
		PORD_INT *	key,
		PORD_INT	scoretype )

Definition at line 536 of file ddcreate.c.

{ PORD_INT *xadj, *adjncy, *vwght, *marker;
  PORD_INT nvtx, nlist, k, weight, deg, u, v, w;
  PORD_INT i, istart, istop, j, jstart, jstop;
 
  nvtx = dd->G->nvtx;
  xadj = dd->G->xadj;
  adjncy = dd->G->adjncy;
  vwght = dd->G->vwght;
  marker = dd->map;
  nlist = nvtx - dd->ndom;
 
  switch(scoretype)
   { case QMRDV:  /* maximal relative decrease of variables in quotient graph */
       for (k = 0; k < nlist; k++)
        { u = msvtxlist[k];
          weight = vwght[u];
          istart = xadj[u];
          istop = xadj[u+1];
          for (i = istart; i < istop; i++)
            weight += vwght[adjncy[i]];
          key[u] = weight / vwght[u];
        }
       break;
 
     case QMD:    /* ----------------------- minimum degree in quotient graph */
       for (k = 0; k < nlist; k++)
        { u = msvtxlist[k];
          marker[u] = -1;
        }
       for (k = 0; k < nlist; k++)
        { u = msvtxlist[k];
          marker[u] = u;
          deg = 0;
          istart = xadj[u];
          istop = xadj[u+1];
          for (i = istart; i < istop; i++)
           { v = adjncy[i];
             jstart = xadj[v];
             jstop = xadj[v+1];
             for (j = jstart; j < jstop; j++)
              { w = adjncy[j];
                if (marker[w] != u)
                 { marker[w] = u;
                   deg += vwght[w];
                 }
              }
           }
          key[u] = deg;
        }
       break;
 
     case QRAND:  /* ------------------------------------------------- random */
       for (k = 0; k < nlist; k++)
        { u = msvtxlist[k];
          key[u] = myrandom(nvtx);
        }
       break;
 
     default:
       fprintf(stderr, "\nError in internal function computePriorities\n"
            "  unrecognized node selection strategy %d\n", scoretype);
       quit();
   }
}

◆ connectedComponents()

PORD_INT connectedComponents ( graph_t * G )

Definition at line 344 of file graph.c.

{ PORD_INT *xadj, *adjncy, *marker, *queue;
  PORD_INT nvtx, u, v, w, qhead, qtail, comp, i, istart, istop;
 
  nvtx = G->nvtx;
  xadj = G->xadj;
  adjncy = G->adjncy;
 
  /* ------------------------
     allocate working storage
     ------------------------ */ 
  mymalloc(marker, nvtx, PORD_INT);
  mymalloc(queue, nvtx, PORD_INT);
 
  /* ---------------
     initializations
     --------------- */
  comp = 0;
  for (u = 0; u < nvtx; u++)
    marker[u] = -1;
 
  /* --------------------------------------
     get the number of connected components
     -------------------------------------- */
  for (u = 0; u < nvtx; u++)
   if (marker[u] == -1)
    { comp++;
      qhead = 0; qtail = 1;
      queue[0] = u; marker[u] = 0;
 
      while (qhead != qtail)       /* breadth first search in each comp. */
       { v = queue[qhead++];
         istart = xadj[v];
         istop = xadj[v+1];
         for (i = istart; i < istop; i++)
          { w = adjncy[i];
            if (marker[w] == -1)
             { queue[qtail++] = w;
               marker[w] = 0;
             }
          }
       }
    }
 
  /* -------------------------------
     free working storage and return
     ------------------------------- */
  free(marker); free(queue);
  return(comp);
}

◆ constructDomainDecomposition()

domdec_t * constructDomainDecomposition	(	graph_t *	G,
		PORD_INT *	map )

Definition at line 472 of file ddcreate.c.

{ domdec_t *dd;
  PORD_INT      *xadj, *adjncy, *vwght, *vtxlist, *vtype, *key, *rep;
  PORD_INT      nvtx, deg, u, i, istart, istop;
 
  nvtx = G->nvtx;
  xadj = G->xadj;
  adjncy = G->adjncy;
  vwght = G->vwght;
 
  /* ---------------------------------------------------------
     sort the vertices in G in ascending order of their degree
     --------------------------------------------------------- */
  mymalloc(vtxlist, nvtx, PORD_INT);
  mymalloc(key, nvtx, PORD_INT);
  for (u = 0; u < nvtx; u++)
   { vtxlist[u] = u;
     istart = xadj[u];
     istop = xadj[u+1];
     switch(G->type)
      { case UNWEIGHTED:
          deg = istop - istart;
          break;
        case WEIGHTED:
          deg = 0;
          for (i = istart; i < istop; i++)
            deg += vwght[adjncy[i]];
          break;
        default:
          fprintf(stderr, "\nError in function constructDomainDecomposition\n"
               "  unrecognized graph type %d\n", G->type);
          quit();
      }
     key[u] = deg;
   }
  distributionCounting(nvtx, vtxlist, key);
  free(key);
 
   /* -------------------------------------------------------------
      build initial domains and cluster multisecs that do not share
      a common domain
      ------------------------------------------------------------- */
   mymalloc(vtype, nvtx, PORD_INT);
   mymalloc(rep, nvtx, PORD_INT);
   for (u = 0; u < nvtx; u++)
    { vtype[u] = 0;
      rep[u] = u;
    }
   buildInitialDomains(G, vtxlist, vtype, rep);
   mergeMultisecs(G, vtype, rep);
   free(vtxlist);
 
   /* --------------------------------------------------
      finally, build the domain decomposition and return
      -------------------------------------------------- */
   dd = initialDomainDecomposition(G, map, vtype, rep);
   free(vtype); free(rep);
   return(dd);
}

◆ constructLevelSep()

void constructLevelSep	(	domdec_t *	dd,
		PORD_INT	domain )

Definition at line 208 of file ddbisect.c.

{ PORD_INT *xadj, *adjncy, *vwght, *vtype, *color, *cwght;
  PORD_INT *queue, *deltaS, *deltaB, *deltaW;
  PORD_INT nvtx, bestvalue, weight, qhead, qtail, qopt, q, dS, dB, dW;
  PORD_INT u, v, w, i, istart, istop, j, jstart, jstop;
 
  /* ======================================================================
     vtype[u]: (u domain)
         1 => domain u has not been touched yet (not in queue, no color flip)
        -1 => domain u is in queue and its deltaS, deltaB, deltaW values
              have to be updated
        -2 => domain u is in queue and no update necessary
        -3 => domain u has flipped its color to black
     deltaS[u], deltaB[u], deltaW[u]:
        u domain: denotes the change in partition size, if u flips its color
        u multisec: deltaB/deltaW denote number of adj. black/white domains
     ====================================================================== */
 
  nvtx = dd->G->nvtx;
  xadj = dd->G->xadj;
  adjncy = dd->G->adjncy;
  vwght = dd->G->vwght;
  vtype = dd->vtype;
  color = dd->color;
  cwght = dd->cwght;
 
  /* ------------------------------------------
     allocate working storage + initializations
     ------------------------------------------ */
  mymalloc(queue, nvtx, PORD_INT);
  mymalloc(deltaS, nvtx, PORD_INT);
  mymalloc(deltaB, nvtx, PORD_INT);
  mymalloc(deltaW, nvtx, PORD_INT);
  for (u = 0; u < nvtx; u++)
   { deltaS[u] = deltaB[u] = deltaW[u] = 0;
     if (vtype[u] == 2)
       deltaW[u] = xadj[u+1] - xadj[u];
   }
 
  /* ---------------------------------------------
     build a BFS tree rooted at domain
     the separator is given by the level structure
     --------------------------------------------- */
  queue[0] = domain;
  qhead = 0; qtail = 1;
  vtype[domain] = -1;
  while ((cwght[BLACK] < cwght[WHITE]) && (qhead != qtail))
   { qopt = 0;
     bestvalue = MAX_INT;
     
     /* --------------------------------------------------------------------
        run through queue, update domains if necessary, and find best domain
        -------------------------------------------------------------------- */
     for (q = qhead; q < qtail; q++)
      { u = queue[q];
        if (vtype[u] == -1)
         { dB = vwght[u]; dW = -dB; dS = 0;
           istart = xadj[u];
           istop = xadj[u+1];
           for (i = istart; i < istop; i++)
            { v = adjncy[i];                    /* color of multisec v */
              weight = vwght[v];                /* is GRAY or WHITE */
              if (color[v] == WHITE)
               { dW -= weight; dS += weight; }  /* multisec will move to S */
              else if (deltaW[v] == 1)
               { dB += weight; dS -= weight; }  /* multisec will move to B */
            }
           deltaS[u] = dS; deltaB[u] = dB; deltaW[u] = dW;
           vtype[u] = -2;
         }
        if (cwght[GRAY] + deltaS[u] < bestvalue)
         { bestvalue = cwght[GRAY] + deltaS[u];
           qopt = q;
         }
      }
 
     /* ----------------------------------------------------
        move best domain to head of queue and color it black
        ---------------------------------------------------- */
     u = queue[qopt];
     swap(queue[qopt], queue[qhead], v);
     qhead++;
     color[u] = BLACK;
     cwght[GRAY] += deltaS[u];
     cwght[BLACK] += deltaB[u];
     cwght[WHITE] += deltaW[u];
     vtype[u] = -3;
 
     /* ------------------------------------------------------------
        update all multisecs that are adjacent to domain u and check
        domains adjacent to the multisecs
        ------------------------------------------------------------ */
     istart = xadj[u];
     istop = xadj[u+1];
     for (i = istart; i < istop; i++)
      { v = adjncy[i];
        deltaB[v]++; deltaW[v]--;
        if (deltaW[v] == 0)       /* color of multisec v changed to BLACK */
          color[v] = BLACK;
        else if (deltaB[v] == 1)  /* color of multisec v changed to GRAY */
         { color[v] = GRAY;
           jstart = xadj[v];
           jstop = xadj[v+1];
           for (j = jstart; j < jstop; j++)
            { w = adjncy[j];
              if (vtype[w] == 1)           /* a new domain enters the queue */
               { queue[qtail++] = w;
                 vtype[w] = -1;
               }
              else if (vtype[w] == -2)     /* update (old) domain in queue */
                vtype[w] = -1;
            }
         }
        else if (deltaW[v] == 1)  /* color of multisec v remains GRAY for */
         { jstart = xadj[v];      /* the last time */
           jstop = xadj[v+1];
           for (j = jstart; j < jstop; j++)
            { w = adjncy[j];
              if (vtype[w] == -2)
                vtype[w] = -1;
            }
         }
      }
   }
 
  /* ---------------------------
     reset vtype and free memory
     --------------------------- */
  for (i = 0; i < qtail; i++)
   { u = queue[i];
     vtype[u] = 1;
   }
  free(queue);
  free(deltaS); free(deltaB); free(deltaW);
}

◆ constructMultisector()

multisector_t * constructMultisector	(	graph_t *	G,
		options_t *	options,
		timings_t *	cpus )

Definition at line 124 of file multisector.c.

{ multisector_t *ms;
  nestdiss_t    *ndroot;
  PORD_INT           *map, nvtx, ordtype;
 
  nvtx = G->nvtx;
 
  /* ------------------------------
     check number of nodes in graph
     ------------------------------ */
  /* -----------------------------------
     JY: inserted the condition
     "&& (options[OPTION_MSGLVL] > 0)"
     below, to avoid systematic printing
     ----------------------------------- */
    if ((nvtx <= MIN_NODES) && (options[OPTION_ORDTYPE] != MINIMUM_PRIORITY)
        && (options[OPTION_MSGLVL] > 0))
   { printf("\nWarning in constructMultisector\n"
         "  graph has less than %d nodes, skipping separator construction\n\n",
         MIN_NODES);
     options[OPTION_ORDTYPE] = MINIMUM_PRIORITY;
   }
  /* --------------------------------------------------------
     determine the multisector according to the ordering type
     -------------------------------------------------------- */
  ordtype = options[OPTION_ORDTYPE];
   switch(ordtype)
   { case MINIMUM_PRIORITY:
       ms = trivialMultisector(G);
       break;
 
     case INCOMPLETE_ND:
     case MULTISECTION:
     case TRISTAGE_MULTISECTION:
       mymalloc(map, nvtx, PORD_INT);
       ndroot = setupNDroot(G, map);
       buildNDtree(ndroot, options, cpus);
       if (ordtype == MULTISECTION)
         ms = extractMS2stage(ndroot);
       else
         ms = extractMSmultistage(ndroot);
       freeNDtree(ndroot);
       freeNDnode(ndroot);
       free(map);
       break;
 
     default:
       fprintf(stderr, "\nError in function constructMultisector\n"
            "  unrecognized ordering type %d\n", ordtype);
       quit();
   }
  return(ms);
}

◆ constructSeparator()

void constructSeparator	(	gbisect_t *	Gbisect,
		options_t *	options,
		timings_t *	cpus )

Definition at line 182 of file gbisect.c.

{ domdec_t *dd, *dd2;
  PORD_INT      *color, *cwght, *map, nvtx, u, i;
 
  nvtx = Gbisect->G->nvtx;
  color = Gbisect->color;
  cwght = Gbisect->cwght;
 
  /* --------------------------------------------------------------
     map vector identifies vertices of Gbisect->G in domain decomp.
     -------------------------------------------------------------- */
  mymalloc(map, nvtx, PORD_INT);
 
  /* --------------------------------------
     construct initial domain decomposition
     -------------------------------------- */
  pord_starttimer(cpus[TIME_INITDOMDEC]);
  dd = constructDomainDecomposition(Gbisect->G, map);
 
#ifdef BE_CAUTIOUS
  checkDomainDecomposition(dd);
#endif
 
  if (options[OPTION_MSGLVL] > 2)
    printf("\t  0. dom.dec.: #nodes %d (#domains %d, weight %d), #edges %d\n",
           dd->G->nvtx, dd->ndom, dd->domwght, dd->G->nedges >> 1);
  pord_stoptimer(cpus[TIME_INITDOMDEC]);
 
  /* ---------------------------------------------------
     construct sequence of coarser domain decompositions
     --------------------------------------------------- */
  pord_starttimer(cpus[TIME_COARSEDOMDEC]);
  i = 0;
  while ((dd->ndom > MIN_DOMAINS) && (i < MAX_COARSENING_STEPS)
         && ((dd->G->nedges >> 1) > dd->G->nvtx))
   { shrinkDomainDecomposition(dd, options[OPTION_NODE_SELECTION3]);
     dd = dd->next; i++;
 
#ifdef BE_CAUTIOUS
     checkDomainDecomposition(dd);
#endif
 
     if (options[OPTION_MSGLVL] > 2)
       printf("\t %2d. dom.dec.: #nodes %d (#domains %d, weight %d), #edges %d"
              "\n", i, dd->G->nvtx, dd->ndom, dd->domwght, dd->G->nedges >> 1);
   }
  pord_stoptimer(cpus[TIME_COARSEDOMDEC]);
 
  /* -----------------------------------------------
     determine coloring of last domain decomposition
     ------------------------------------------------ */
  pord_starttimer(cpus[TIME_INITSEP]);
  initialDDSep(dd);
  if (dd->cwght[GRAY] > 0)
    improveDDSep(dd);
 
#ifdef BE_CAUTIOUS
  checkDDSep(dd);
#endif
 
  if (options[OPTION_MSGLVL] > 2)
    printf("\t %2d. dom.dec. sep.: S %d, B %d, W %d [cost %7.2f]\n",
           i, dd->cwght[GRAY], dd->cwght[BLACK], dd->cwght[WHITE],
           F(dd->cwght[GRAY], dd->cwght[BLACK], dd->cwght[WHITE]));
  pord_stoptimer(cpus[TIME_INITSEP]);
 
  /* --------------
     refine coloring
     --------------- */
 
  pord_starttimer(cpus[TIME_REFINESEP]);
  while (dd->prev != NULL)
   { dd2 = dd->prev;
     dd2->cwght[GRAY] = dd->cwght[GRAY];
     dd2->cwght[BLACK] = dd->cwght[BLACK];
     dd2->cwght[WHITE] = dd->cwght[WHITE];
     for (u = 0; u < dd2->G->nvtx; u++)
       dd2->color[u] = dd->color[dd2->map[u]];
     freeDomainDecomposition(dd);
     if (dd2->cwght[GRAY] > 0)
       improveDDSep(dd2);
 
#ifdef BE_CAUTIOUS
     checkDDSep(dd2);
#endif
 
     dd = dd2;
     i--;
     if (options[OPTION_MSGLVL] > 2)
       printf("\t %2d. dom.dec. sep.: S %d, B %d, W %d [cost %7.2f]\n",
              i, dd->cwght[GRAY], dd->cwght[BLACK], dd->cwght[WHITE],
              F(dd->cwght[GRAY], dd->cwght[BLACK], dd->cwght[WHITE]));
   }
  pord_stoptimer(cpus[TIME_REFINESEP]);
 
  /* ---------------------------------
     copy coloring to subgraph Gbisect
     --------------------------------- */
  cwght[GRAY] = dd->cwght[GRAY];
  cwght[BLACK] = dd->cwght[BLACK];
  cwght[WHITE] = dd->cwght[WHITE];
  for (u = 0; u < nvtx; u++)
    color[u] = dd->color[map[u]];
  freeDomainDecomposition(dd);
  free(map);
}

◆ crunchElimGraph()

PORD_INT crunchElimGraph ( gelim_t * Gelim )

Definition at line 302 of file gelim.c.

{ PORD_INT *xadj, *adjncy, *len;
  PORD_INT nvtx, nedges, u, i, isrc, idest;
 
  nvtx = Gelim->G->nvtx;
  nedges = Gelim->G->nedges;
  xadj = Gelim->G->xadj;
  adjncy = Gelim->G->adjncy;
  len = Gelim->len;
 
  /* ---------------------------------------------
     mark begining of u's adjacency list by -(u+1)
     --------------------------------------------- */
  for (u = 0; u < nvtx; u++)
   { i = xadj[u];                /* is adjacency list of u still in use? */
     if (i != -1)                /* verify that list is non-empty */
      { if (len[u] == 0)
         { fprintf(stderr, "\nError in function crunchElimGraph\n"
                "  adjacency list of node %d is empty\n", u);
           quit();
         }
        xadj[u] = adjncy[i];     /* if so, move first item to xadj[u] */
        adjncy[i] = -(u+1);      /* u's adjacency list is headed by -(u+1) */
        if (len[u] == 0)
          printf("error: u %d, len %d\n", u, len[u]);
      }
   }
 
  /* --------------------------
     crunch all adjacency lists
     -------------------------- */
  idest = isrc = 0;
  while (isrc < Gelim->G->nedges)
   { u = adjncy[isrc++];
     if (u < 0)                        /* a new adjacency list starts here */
      { u = -u - 1;                    /* it's the adjacency list of u */
        adjncy[idest] = xadj[u];       /* first item was stored in xadj[u] */
        xadj[u] = idest++;
        for (i = 1; i < len[u]; i++)
          adjncy[idest++] = adjncy[isrc++];
      }
   }
  Gelim->G->nedges = idest;
 
  /* ------------------
     was it successful?
     ------------------ */
  if (idest < nedges) return(TRUE);
  else return (FALSE);
}

◆ distributionCounting()

void distributionCounting	(	PORD_INT	n,
		PORD_INT *	node,
		PORD_INT *	key )

Definition at line 186 of file sort.c.

{ register PORD_INT i;
  PORD_INT *tmp, *count, minkey, maxkey, l, u, vk;
 
  /* determine maximal and minimal key */
  minkey = MAX_INT;
  maxkey = 0;
  for (i = 0; i < n; i++)
   { u = node[i];
     maxkey = max(key[u], maxkey);
     minkey = min(key[u], minkey);
   }
  l = maxkey-minkey;
  /* printf("minkey %d, maxkey %d, range %d\n", minkey, maxkey, l); */
  mymalloc(count, (l+1), PORD_INT);
  mymalloc(tmp, n, PORD_INT);
  for (i = 0; i <= l; i++)
    count[i] = 0;
 
  /* scale down all key-values */
  for (i = 0; i < n; i++)
   { u = node[i];
     vk = key[u]-minkey;
     key[u] = vk;
     count[vk]++;
   }
 
  /* now do the sorting */
  for (i = 1; i <= l; i++)
    count[i] += count[i-1];
  for (i = n-1; i >= 0; i--)
   { u = node[i];
     tmp[--count[key[u]]] = u;
   }
  for (i = 0; i < n; i++)
    node[i] = tmp[i];
/*
  for (i = 0; i < n; i++)
   { u = node[i];
     printf("  node %d, key %d\n", u, key[u]);
   }
*/
  free(count);
  free(tmp);
} 

◆ DMviaFlow()

void DMviaFlow	(	gbipart_t *	Gbipart,
		PORD_INT *	flow,
		PORD_INT *	rc,
		PORD_INT *	dmflag,
		PORD_INT *	dmwght )

Definition at line 528 of file gbipart.c.

{ PORD_INT *xadj, *adjncy, *vwght, *queue, qhead, qtail;
  PORD_INT u, v, x, nX, y, nY, i, istart, istop;
 
  xadj = Gbipart->G->xadj;
  adjncy = Gbipart->G->adjncy;
  vwght = Gbipart->G->vwght;
  nX = Gbipart->nX;
  nY = Gbipart->nY;
 
  mymalloc(queue, (nX+nY), PORD_INT);
 
  /* ----------------------------------------------------------
     mark all nodes reachable from source/sink with SOURCE/SINK
     ---------------------------------------------------------- */
  qhead = qtail = 0;
  for (x = 0; x < nX; x++)
    if (rc[x] > 0)
     { queue[qtail++] = x;
       dmflag[x] = SOURCE;
     }
    else dmflag[x] = FREE;
  for (y = nX; y < nX+nY; y++)
    if (rc[y] > 0)
     { queue[qtail++] = y;
       dmflag[y] = SINK;
     }
    else dmflag[y] = FREE;
 
  /* --------------------------------------------------------------------
     construct Dulmage-Mendelsohn decomp. starting with SOURCE/SINK nodes
     -------------------------------------------------------------------- */
  while (qhead != qtail)
   { u = queue[qhead++];
     istart = xadj[u];
     istop = xadj[u+1];
     switch(dmflag[u])
      { case SOURCE:
          for (i = istart; i < istop; i++)
           { v = adjncy[i];
             if ((dmflag[v] == FREE) && ((v >= nX) || (flow[i] < 0)))
              { queue[qtail++] = v;
                dmflag[v] = SOURCE;    /* v reachable via forward edge u->v */
              }                        /*         or via backward edge u<-v */
           }
          break;
        case SINK:
          for (i = istart; i < istop; i++)
           { v = adjncy[i];
             if ((dmflag[v] == FREE) && ((v < nX) || (flow[i] > 0)))
              { queue[qtail++] = v;
                dmflag[v] = SINK;    /* u reachable via forward edge v->u */
              }                      /*         or via backward edge v<-u */
           }
          break;
      }
   }
 
  /* -----------------------------------------------------
     all nodes x in X with dmflag[x] = SOURCE belong to SI
     all nodes x in X with dmflag[x] = SINK   belong to SX
     all nodes x in X with dmflag[x] = FREE   belong to SR
     ----------------------------------------------------- */
  dmwght[SI] = dmwght[SX] = dmwght[SR] = 0;
  for (x = 0; x < nX; x++)
    switch(dmflag[x])
     { case SOURCE: dmflag[x] = SI; dmwght[SI] += vwght[x]; break;
       case SINK:   dmflag[x] = SX; dmwght[SX] += vwght[x]; break;
       default:     dmflag[x] = SR; dmwght[SR] += vwght[x];
     }
 
  /* -----------------------------------------------------
     all nodes y in Y with dmflag[y] = SOURCE belong to BX
     all nodes y in Y with dmflag[y] = SINK   belong to BI
     all nodes y in Y with dmflag[y] = FREE   belong to BR
     ----------------------------------------------------- */
  dmwght[BI] = dmwght[BX] = dmwght[BR] = 0;
  for (y = nX; y < nX+nY; y++)
    switch(dmflag[y])
     { case SOURCE: dmflag[y] = BX; dmwght[BX] += vwght[y]; break;
       case SINK:   dmflag[y] = BI; dmwght[BI] += vwght[y]; break;
       default:     dmflag[y] = BR; dmwght[BR] += vwght[y];
     }
 
  free(queue);
}

◆ DMviaMatching()

void DMviaMatching	(	gbipart_t *	Gbipart,
		PORD_INT *	matching,
		PORD_INT *	dmflag,
		PORD_INT *	dmwght )

Definition at line 435 of file gbipart.c.

{ PORD_INT *xadj, *adjncy, *vwght, *queue, qhead, qtail;
  PORD_INT u, x, nX, y, nY, i, istart, istop;
 
  xadj = Gbipart->G->xadj;
  adjncy = Gbipart->G->adjncy;
  vwght = Gbipart->G->vwght;
  nX = Gbipart->nX;
  nY = Gbipart->nY;
 
  mymalloc(queue, (nX+nY), PORD_INT);
 
  /* ----------------------------------------------------------------------
     mark all exposed nodes of X with SI and all exposed nodes of Y with BI
     ---------------------------------------------------------------------- */
  qhead = qtail = 0;
  for (x = 0; x < nX; x++)
    if (matching[x] == FREE)
     { queue[qtail++] = x;
       dmflag[x] = SI;
     }
    else dmflag[x] = SR;
  for (y = nX; y < nX+nY; y++)
    if (matching[y] == FREE)
     { queue[qtail++] = y;
       dmflag[y] = BI;
     }
    else dmflag[y] = BR;
 
  /* ------------------------------------------------------------------
     construct Dulmage-Mendelsohn decomp. starting with SI and BI nodes
     ------------------------------------------------------------------ */
  while (qhead != qtail)
   { u = queue[qhead++];
     istart = xadj[u];
     istop = xadj[u+1];
     switch(dmflag[u])
      { case SI:
          for (i = istart; i < istop; i++)
           { y = adjncy[i];
             if (dmflag[y] == BR)
              { queue[qtail++] = y;
                dmflag[y] = BX;
              }
           }
          break;
        case BX:
          x = matching[u];
          dmflag[x] = SI;
          queue[qtail++] = x;
          break;
        case BI:
          for (i = istart; i < istop; i++)
           { x = adjncy[i];
             if (dmflag[x] == SR)
              { queue[qtail++] = x;
                dmflag[x] = SX;
              }
           }
          break;
        case SX:
          y = matching[u];
          dmflag[y] = BI;
          queue[qtail++] = y;
          break;
      }
   }
 
  /* ----------------------
     fill the dmwght vector
     ---------------------- */
  dmwght[SI] = dmwght[SX] = dmwght[SR] = 0;
  for (x = 0; x < nX; x++)
    switch(dmflag[x])
     { case SI: dmwght[SI] += vwght[x]; break;
       case SX: dmwght[SX] += vwght[x]; break;
       case SR: dmwght[SR] += vwght[x]; break;
     }
  dmwght[BI] = dmwght[BX] = dmwght[BR] = 0;
  for (y = nX; y < nX+nY; y++)
    switch(dmflag[y])
     { case BI: dmwght[BI] += vwght[y]; break;
       case BX: dmwght[BX] += vwght[y]; break;
       case BR: dmwght[BR] += vwght[y]; break;
     }
 
  free(queue);
}

◆ eliminateMultisecs()

void eliminateMultisecs	(	domdec_t *	dd,
		PORD_INT *	msvtxlist,
		PORD_INT *	rep )

Definition at line 606 of file ddcreate.c.

{ PORD_INT *xadj, *adjncy, *vtype;
  PORD_INT nvtx, nlist, keepon, u, v, w, k, i, istart, istop;
 
  nvtx = dd->G->nvtx;
  xadj = dd->G->xadj;
  adjncy = dd->G->adjncy;
  vtype = dd->vtype;
  nlist = nvtx - dd->ndom;
 
  /* -------------------------------------------------------
     eliminate multisecs according to the order in msvtxlist
     ------------------------------------------------------- */
  for (k = 0; k < nlist; k++)
   { u = msvtxlist[k];
     istart = xadj[u];
     istop = xadj[u+1];
     keepon = TRUE;
     for (i = istart; i < istop; i++)
      { v = adjncy[i];
        if (rep[v] != v)        /* domain already absorbed by an eliminated */
         { keepon = FALSE;      /* multisec => multisec u cannot be deleted */
           break;
         }
      }
     if (keepon)
      { vtype[u] = 3;
        for (i = istart; i < istop; i++)
         { v = adjncy[i];
           rep[v] = u;
         }
      }
   }
 
  /* ------------------------------------------------------------
     eliminate all multisecs that are adjacent to only one domain
     ------------------------------------------------------------ */
  for (k = 0; k < nlist; k++)
   { u = msvtxlist[k];
     if (vtype[u] == 2)
      { v = -1;
        istart = xadj[u];
        istop = xadj[u+1];
        for (i = istart; i < istop; i++)
         { w = adjncy[i];
           if (v == -1)
             v = rep[w];          /* u adjacent to domain v = rep[w] */
           else if (v != rep[w])
            { v = -1;             /* u adjacent to another domain */
              break;
            }
         }
        if (v != -1)              /* u absorbed by domain v */
         { vtype[u] = 4;
           rep[u] = v;
         }
      }
   }
}

◆ eliminateStage()

void eliminateStage	(	minprior_t *	minprior,
		PORD_INT	istage,
		PORD_INT	scoretype,
		timings_t *	cpus )

Definition at line 232 of file minpriority.c.

{ gelim_t     *Gelim;
  bucket_t    *bucket;
  stageinfo_t *stageinfo;
  PORD_INT         *stage, *reachset, *auxbin, *auxtmp, *auxaux;
  PORD_INT         *degree, *score;
  PORD_INT         *pflag, nreach, nvtx, r, u, i;
 
  Gelim = minprior->Gelim;
  bucket = minprior->bucket;
  stage = minprior->ms->stage;
  stageinfo = minprior->stageinfo + istage;
  reachset = minprior->reachset;
  auxaux = minprior->auxaux;
  auxbin = minprior->auxbin;
  auxtmp = minprior->auxtmp;
  pflag = &(minprior->flag);
  
  nvtx = Gelim->G->nvtx;
  degree = Gelim->degree;
  score = Gelim->score;
 
#ifdef DEBUG
  printf("\nSTARTING NEW ELIMINATION STAGE (nedges %d, maxedges %d)\n\n",
         Gelim->G->nedges, Gelim->maxedges);
  if (istage> 0) printElimGraph(Gelim);
  /* waitkey(); */
#endif
 
  /* -------------------------------------------------------------
     load reachset with all principal variables in stage <= istage
     ------------------------------------------------------------- */
  nreach = 0;
  for (u = 0; u < nvtx; u++)
    if ((score[u] == -1) && (stage[u] <= istage))
     { reachset[nreach++] = u;
       score[u] = degree[u];
       /* score[u] = degree[u]*(degree[u]-1)/2; */
     }
 
  /* ----------------------------------------------------------------
     do an initial update of the vertices in reachset and fill bucket
     ---------------------------------------------------------------- */
  pord_starttimer(cpus[TIME_UPDSCORE]);
  updateDegree(Gelim, reachset, nreach, auxbin);
  updateScore(Gelim, reachset, nreach, scoretype, auxbin);
  pord_stoptimer(cpus[TIME_UPDSCORE]);
  for (i = 0; i < nreach; i++)
   { u = reachset[i];
     insertBucket(bucket, score[u], u);
   }
 
  /* -------------------------------------
     and now start the elimination process
     ------------------------------------- */
  while (TRUE)
   { if (eliminateStep(minprior, istage, scoretype) == 0)
       break;
     nreach = minprior->nreach;
 
#ifdef BE_CAUTIOUS
     printf("checking arrays auxtmp and auxbin\n");
     for (u = 0; u < nvtx; u++)
       if ((auxtmp[u] >= *pflag) || (auxbin[u] != -1))
        { printf("ERROR: flag = %d, auxtmp[%d] = %d, auxbin[%d] = %d\n",
                 *pflag, u, auxtmp[u], u, auxbin[u]);
          quit();
        }
#endif
 
     /* ----------------------------------------------------------
        update the adjacency structure of all vertices in reachset
        ---------------------------------------------------------- */
     pord_starttimer(cpus[TIME_UPDADJNCY]);
     updateAdjncy(Gelim, reachset, nreach, auxtmp, pflag);
     pord_stoptimer(cpus[TIME_UPDADJNCY]);
 
     /* ----------------------------------------
        find indistinguishable nodes in reachset
        ---------------------------------------- */
     pord_starttimer(cpus[TIME_FINDINODES]);
     findIndNodes(Gelim, reachset, nreach, auxbin, auxaux, auxtmp, pflag);
     pord_stoptimer(cpus[TIME_FINDINODES]);
 
#ifdef BE_CAUTIOUS
     printf("checking arrays auxtmp and auxbin\n");
     for (u = 0; u < nvtx; u++)
       if ((auxtmp[u] >= *pflag) || (auxbin[u] != -1))
        { printf("ERROR: flag = %d, auxtmp[%d] = %d, auxbin[%d] = %d\n",
                 *pflag, u, auxtmp[u], u, auxbin[u]);
          quit();
        }
#endif
 
     /* ----------------------------------------------------------------
        clean reachset of nonprincipal nodes and nodes not in this stage
        ---------------------------------------------------------------- */
     r = 0;
     for (i = 0; i < nreach; i++)
      { u = reachset[i];
        if (score[u] >= 0)
          reachset[r++] = u;
      }
     nreach = r;
 
     /* ---------------------------------------------------
        update the degree/score of all vertices in reachset
        --------------------------------------------------- */
     pord_starttimer(cpus[TIME_UPDSCORE]);
     updateDegree(Gelim, reachset, nreach, auxbin);
     updateScore(Gelim, reachset, nreach, scoretype, auxbin);
     pord_stoptimer(cpus[TIME_UPDSCORE]);
 
     /* ----------------------------
        re-insert vertices in bucket
        ---------------------------- */
     for (i = 0; i < nreach; i++)
      { u = reachset[i];
        insertBucket(bucket, score[u], u);
      }
 
     stageinfo->nstep++;
   }
}

◆ eliminateStep()

PORD_INT eliminateStep	(	minprior_t *	minprior,
		PORD_INT	istage,
		PORD_INT	scoretype )

Definition at line 361 of file minpriority.c.

{ gelim_t     *Gelim;
  bucket_t    *bucket;
  stageinfo_t *stageinfo;
  PORD_INT         *stage, *reachset, *auxtmp;
  PORD_INT         *xadj, *adjncy, *vwght, *len, *degree, *score;
  PORD_INT         *pflag, *pnreach, nelim, minscr, vwghtu, u, v, i, istart, istop;
  FLOAT       tri, rec;
 
  Gelim = minprior->Gelim;
  bucket = minprior->bucket;
  stage = minprior->ms->stage;
  stageinfo = minprior->stageinfo + istage;
  reachset = minprior->reachset;
  pnreach = &(minprior->nreach);
  auxtmp = minprior->auxtmp;
  pflag = &(minprior->flag);
 
  xadj = Gelim->G->xadj;
  adjncy = Gelim->G->adjncy;
  vwght = Gelim->G->vwght;
  len = Gelim->len;
  degree = Gelim->degree;
  score = Gelim->score;
 
#ifdef DEBUG
  printf("\nStarting new elimination step (nedges %d, maxedges %d)\n",
         Gelim->G->nedges, Gelim->maxedges);
  /* waitkey(); */
#endif
 
  /* ----------------------
     check for empty bucket
     ---------------------- */
  if ((u = minBucket(bucket)) == -1)
    return(0);
  minscr = score[u];
 
  /* ----------------------------------------
     loop while nodes of minimum score remain
     ---------------------------------------- */
  nelim = 0;
  *pnreach = 0;
  while (TRUE)
   { vwghtu = vwght[u];
 
     /* --------------------------------------------------
        increment welim and nelim and remove u from bucket
        -------------------------------------------------- */
     removeBucket(bucket, u);
     stageinfo->welim += vwghtu;
     nelim++;
 
     /* -----------------------------------------------------------------
        call buildElement to create element u and merge u's boundary with 
        the nodes in reachset; remove any vertex from bucket that belongs
        to u's boundary and to the actual stage
        ----------------------------------------------------------------- */
     buildElement(Gelim, u);
     istart = xadj[u];
     istop = istart + len[u];
     for (i = istart; i < istop; i++)
      { v = adjncy[i];                 /* v belongs to u's boundary */
        if (auxtmp[v] < *pflag)        /* v not yet in reachset */
         { auxtmp[v] = *pflag;
           if (stage[v] <= istage)     /* v belongs to actual stage */
             removeBucket(bucket, v);
           reachset[(*pnreach)++] = v;
         }
      }
 
#ifdef DEBUG
     printf("Node %d (weight %d, score %d) eliminated: (boundary weight %d)\n",
            u, vwghtu, minscr, degree[u]);
     for (i = istart; i < istop; i++)
       printf("%4d (degree %2d)", adjncy[i], degree[adjncy[i]]);
     printf("\n");
#endif
 
      /* ---------------------------------------------------------------
         increment the storage and operation counts for this elim. stage
         --------------------------------------------------------------- */
     tri = vwghtu;
     rec = degree[u];
     stageinfo->nzf += (PORD_INT)((tri * (tri+1)) / 2);
     stageinfo->nzf += (PORD_INT)(tri * rec);
     stageinfo->ops += (tri*tri*tri) / 3.0 + (tri*tri) / 2.0 - (5*tri) / 6.0;
     stageinfo->ops += (tri*tri*rec) + (rec*(rec+1)*tri);
 
     /* ---------------------------------------------------------------
        end this elim. step, if one of the following conditions is true
         (1) no multiple elimination
         (2) bucket empty
         (3) no further variable with minimum score
        ---------------------------------------------------------------- */
     if (scoretype / 10 == 0)
       break;
     if ((u = minBucket(bucket)) == -1)
       break;
     if (score[u] > minscr)
       break;
   }
 
  /* -----------------------
     clear auxtmp and return
     ----------------------- */
  (*pflag)++;
  return(nelim);
}

◆ exchangeBuffer()

buffer_t * exchangeBuffer	(	topology_t *	,
		buffer_t *	,
		PORD_INT	)

◆ expandElimTree()

elimtree_t * expandElimTree	(	elimtree_t *	T,
		PORD_INT *	vtxmap,
		PORD_INT	nvtxorg )

Definition at line 517 of file tree.c.

{ elimtree_t *T2;
  PORD_INT        *vtx2front, *vtx2front2;
  PORD_INT        nfronts, J, u;
 
  nfronts = T->nfronts;
 
  /* --------------------------------------------------------------
     allocate space for the new elimtree object and copy front data
     the expanded tree has the same number of fronts/tree structure
     -------------------------------------------------------------- */
  T2 = newElimTree(nvtxorg, nfronts);
  T2->root = T->root;
  for (J = 0; J < nfronts; J++)
   { T2->ncolfactor[J] = T->ncolfactor[J];
     T2->ncolupdate[J] = T->ncolupdate[J];
     T2->parent[J] = T->parent[J];
     T2->firstchild[J] = T->firstchild[J];
     T2->silbings[J] = T->silbings[J];
   }
 
  /* ---------------------------------------------------------------------
     set up the new vtx2front vector; the trees only differ in this vector
     --------------------------------------------------------------------- */
  vtx2front = T->vtx2front;
  vtx2front2 = T2->vtx2front;
  for (u = 0; u < nvtxorg; u++)
    vtx2front2[u] = vtx2front[vtxmap[u]];
 
  return(T2);
}

◆ extendedAdd()

denseMtx_t * extendedAdd	(	denseMtx_t *	,
		denseMtx_t *	,
		PORD_INT *	,
		PORD_INT *	)

◆ extendedAddPAR()

denseMtx_t * extendedAddPAR	(	denseMtx_t *	,
		denseMtx_t *	,
		PORD_INT *	,
		PORD_INT *	)

◆ extractElimTree()

elimtree_t * extractElimTree ( gelim_t * Gelim )

Definition at line 1012 of file gelim.c.

{ elimtree_t *T;
  PORD_INT        *vwght, *par, *degree, *score, *sib, *fch;
  PORD_INT        *ncolfactor, *ncolupdate, *parent, *vtx2front;
  PORD_INT        nvtx, nfronts, root, u, v, front;
 
  nvtx = Gelim->G->nvtx;
  vwght = Gelim->G->vwght;
  par = Gelim->parent;
  degree = Gelim->degree;
  score = Gelim->score;
 
  /* ------------------------
     allocate working storage
     ------------------------ */
  mymalloc(sib, nvtx, PORD_INT);
  mymalloc(fch, nvtx, PORD_INT);
  for (u = 0; u < nvtx; u++)
    sib[u] = fch[u] = -1;
 
  /* --------------------------------------------------------------
     count fronts and create top-down view of the tree given by par
     -------------------------------------------------------------- */
  nfronts = 0;
  root = -1;
  for (u = 0; u < nvtx; u++)
    switch(score[u])
     { case -2:          /* variable u is nonprincipal */
         break;
       case -3:          /* variable u has been elim. and now forms an elem. */
         sib[u] = root;
         root = u;
         nfronts++;
         break;
       case -4:          /* element u has been absorbed by par[u] */
         v = par[u];
         sib[u] = fch[v];
         fch[v] = u;
         nfronts++;
         break;
       default:
         fprintf(stderr, "\nError in function extractElimTree\n"
              "  ordering not complete (score[%d] = %d)\n", u, score[u]);
         quit();
     }
 
#ifdef DEBUG
  for (u = 0; u < nvtx; u++)
    printf("node %d: score %d, par %d, fch %d, sib %d\n", u, score[u],
           par[u], fch[u], sib[u]);
#endif
 
  /* --------------------------------------
     allocate space for the elimtree object
     -------------------------------------- */
  T = newElimTree(nvtx, nfronts);
  ncolfactor = T->ncolfactor;
  ncolupdate = T->ncolupdate;
  parent = T->parent;
  vtx2front = T->vtx2front;
 
  /* -------------------------------------------------------------
     fill the vtx2front vector so that representative vertices are
     mapped in a post-order traversal
     ------------------------------------------------------------- */
  nfronts = 0;
  u = root;
  while (u != -1)
   { while (fch[u] != -1)
       u = fch[u];
     vtx2front[u] = nfronts++;
     while ((sib[u] == -1) && (par[u] != -1))
      { u = par[u];
        vtx2front[u] = nfronts++;
      }
     u = sib[u];
   }
 
  /* ---------------------------------------------------
     fill in the vtx2front map for nonprincipal vertices
     --------------------------------------------------- */
  for (u = 0; u < nvtx; u++)
    if (score[u] == -2)
     { v = u;
       while ((par[v] != -1) && (score[v] == -2))
         v = par[v];
       vtx2front[u] = vtx2front[v];
     }
 
  /* -------------------------------------------------------------
     set up the parent vector of T and fill ncolfactor, ncolupdate
     ------------------------------------------------------------- */
  for (u = 0; u < nvtx; u++)
   { front = vtx2front[u];
     if (score[u] == -3)
      { parent[front] = -1;
        ncolfactor[front] = vwght[u];
        ncolupdate[front] = degree[u];
      }
     if (score[u] == -4)
      { parent[front] = vtx2front[par[u]];
        ncolfactor[front] = vwght[u];
        ncolupdate[front] = degree[u];
      }
   }
 
  /* ----------------------------
     set up all other arrays of T
     ---------------------------- */
  initFchSilbRoot(T);
 
    /* ----------------------
     free memory and return
     ---------------------- */
  free(sib); free(fch);
  return(T);
}

◆ extractMS2stage()

multisector_t * extractMS2stage ( nestdiss_t * ndroot )

Definition at line 182 of file multisector.c.

{ multisector_t *ms;
  nestdiss_t    *nd, *parent;
  PORD_INT           *stage, *intvertex, *intcolor;
  PORD_INT           nvint, nnodes, totmswght, i;
 
  /* ----------------------------------------------------------------- 
     allocate memory for the multisector object and init. stage vector
     ----------------------------------------------------------------- */
  ms = trivialMultisector(ndroot->G);
  stage = ms->stage;
 
  /* ------------------------------------------------------------
     extract the stages of the separator vertices:
     stage[u] = 1, iff u belongs to a separator
     ------------------------------------------------------------ */
  nnodes = totmswght = 0;
  for (nd = ndroot; nd->childB != NULL; nd = nd->childB);
  while (nd != ndroot)
   { parent = nd->parent;
     if ((parent == NULL) || (parent->childB == NULL)
        || (parent->childW == NULL))
       { fprintf(stderr, "\nError in function extractMS2stage\n"
              "  nested dissection tree corrupted\n");
         quit();
       }
      if (parent->childB == nd)        /* left subtree of parent visited */
        for (nd = parent->childW; nd->childB != NULL; nd = nd->childB);
      else                             /* right subtree of parent visited */
       { nd = parent;                  /* extract the separator of parent */
         totmswght += nd->cwght[GRAY];
         nvint = nd->nvint;
         intvertex = nd->intvertex;
         intcolor = nd->intcolor;
         for (i = 0; i < nvint; i++)
           if (intcolor[i] == GRAY)
            { nnodes++;
              stage[intvertex[i]] = 1;
            }
       }
   }
 
  /* ------------------------------------------
     finalize the multisector object and return
     ------------------------------------------ */
  ms->nstages = 2;
  ms->nnodes = nnodes;
  ms->totmswght = totmswght;
 
  return(ms);
}

◆ extractMSmultistage()

multisector_t * extractMSmultistage ( nestdiss_t * ndroot )

Definition at line 238 of file multisector.c.

{ multisector_t *ms;
  nestdiss_t    *nd, *parent;
  PORD_INT           *stage, *intvertex, *intcolor;
  PORD_INT           nvtx, nvint, maxstage, istage, nnodes, totmswght, i, u;
 
  /* -----------------------------------------------------------------
     allocate memory for the multisector object and init. stage vector
     ----------------------------------------------------------------- */
  ms = trivialMultisector(ndroot->G);
  stage = ms->stage;
 
  /* ------------------------------------------------------------
     extract the stages of the separator vertices:
     stage[u] = i, i>0, iff u belongs to a separator in depth i-1
     ------------------------------------------------------------ */
  maxstage = nnodes = totmswght = 0;
  for (nd = ndroot; nd->childB != NULL; nd = nd->childB);
  while (nd != ndroot)
   { parent = nd->parent;
     if ((parent == NULL) || (parent->childB == NULL)
        || (parent->childW == NULL))
       { fprintf(stderr, "\nError in function extractMSmultistage\n"
              "  nested dissection tree corrupted\n");
         quit();
       }
      if (parent->childB == nd)        /* left subtree of parent visited */
        for (nd = parent->childW; nd->childB != NULL; nd = nd->childB);
      else                             /* right subtree of parent visited */
       { nd = parent;                  /* extract the separator of parent */
         istage = nd->depth + 1;       /* sep. vertices belong to this stage */
         maxstage = max(maxstage, istage);
         totmswght += nd->cwght[GRAY];
         nvint = nd->nvint;
         intvertex = nd->intvertex;
         intcolor = nd->intcolor;
         for (i = 0; i < nvint; i++)
           if (intcolor[i] == GRAY)
            { nnodes++;
              stage[intvertex[i]] = istage;
            }
       }
   }
 
  /* --------------------------------------------------------------------
     we have: stage[u] = 0 => u belongs to a domain
              stage[u] = 1 => u belongs to the root separator (depth = 0)
                 :
              stage[u] = maxstage => u belongs to a leaf separator
     but we must eliminate the separators in a bottom-up fashion; we like
     to have: stage[u] = 0 => u belongs to a domain
              stage[u] = 1 => u belongs to a leaf separator
                 :
              stage[u] = maxstage => u belongs to the root separator
     -------------------------------------------------------------------- */
  nvtx = ndroot->G->nvtx;
  for (u = 0; u < nvtx; u++)
    if (stage[u] > 0)
      stage[u] = maxstage - stage[u] + 1;
 
  /* ------------------------------------------
     finalize the multisector object and return
     ------------------------------------------ */
  ms->nstages = maxstage + 1;
  ms->nnodes = nnodes;
  ms->totmswght = totmswght;
 
  return(ms);
}

◆ factorize1x1Kernel()

denseMtx_t * factorize1x1Kernel	(	denseMtx_t *	,
		PORD_INT	)

◆ factorize1x1KernelPAR()

denseMtx_t * factorize1x1KernelPAR	(	topology_t *	,
		mask_t *	,
		PORD_INT	,
		denseMtx_t *	,
		frontsub_t *	,
		timings_t *	)

◆ factorize2x2Kernel()

denseMtx_t * factorize2x2Kernel	(	denseMtx_t *	,
		PORD_INT	)

◆ factorize2x2KernelPAR()

denseMtx_t * factorize2x2KernelPAR	(	topology_t *	,
		mask_t *	,
		PORD_INT	,
		denseMtx_t *	,
		frontsub_t *	,
		timings_t *	)

◆ factorize3x3Kernel()

denseMtx_t * factorize3x3Kernel	(	denseMtx_t *	,
		PORD_INT	)

◆ factorize3x3KernelPAR()

denseMtx_t * factorize3x3KernelPAR	(	topology_t *	,
		mask_t *	,
		PORD_INT	,
		denseMtx_t *	,
		frontsub_t *	,
		timings_t *	)

◆ findIndMultisecs()

void findIndMultisecs	(	domdec_t *	dd,
		PORD_INT *	msvtxlist,
		PORD_INT *	rep )

Definition at line 670 of file ddcreate.c.

{ PORD_INT *xadj, *adjncy, *vtype, *tmp, *bin, *checksum, *next, *key;
  PORD_INT nvtx, nlist, flag, keepon, deg, chk, ulast, u, v, k, i, istart, istop;
 
  nvtx = dd->G->nvtx;
  xadj = dd->G->xadj;
  adjncy = dd->G->adjncy;
  vtype = dd->vtype;
  nlist = nvtx - dd->ndom;
  checksum = dd->map;
 
  /* ------------------------
     allocate working storage
     ------------------------ */
  mymalloc(tmp, nvtx, PORD_INT);
  mymalloc(bin, nvtx, PORD_INT);
  mymalloc(next, nvtx, PORD_INT);
  mymalloc(key, nvtx, PORD_INT);
  for (u = 0; u < nvtx; u++)
   { tmp[u] = -1;
     bin[u] = -1;
   }
 
  /* -------------------------------------------------------------------
     compute checksums for all unelim./unabsorbed multisecs in msvtxlist
     ------------------------------------------------------------------- */
  flag = 1;
  for (k = 0; k < nlist; k++)
   { u = msvtxlist[k];
     if (vtype[u] == 2)
      { deg = chk = 0;
        istart = xadj[u];
        istop = xadj[u+1];
        for (i = istart; i < istop; i++)
         { v = adjncy[i];
           if (tmp[rep[v]] != flag)
            { tmp[rep[v]] = flag;
              chk += rep[v];
              deg++;
            }
         }
        chk = chk % nvtx;
        checksum[u] = chk;
        key[u] = deg;
        next[u] = bin[chk];
        bin[chk] = u;
        flag++;
      }
   }
 
  /* ---------------------------------
     merge indistinguishable multisecs
     --------------------------------- */
  for (k = 0; k < nlist; k++)
   { u = msvtxlist[k];
     if (vtype[u] == 2)
      { chk = checksum[u];
        v = bin[chk];          /* examine all multisecs in bin[hash] */
        bin[chk] = -1;         /* do this only once */
        while (v != -1)
         { istart = xadj[v];
           istop = xadj[v+1];
           for (i = istart; i < istop; i++)
             tmp[rep[adjncy[i]]] = flag;
           ulast = v;          /* v is principal and u is a potiential */
           u = next[v];        /* nonprincipal variable */
           while (u != -1)
            { keepon = TRUE;
              if (key[u] != key[v])
                keepon = FALSE;
              if (keepon)
               { istart = xadj[u];
                 istop = xadj[u+1];
                 for (i = istart; i < istop; i++)
                   if (tmp[rep[adjncy[i]]] != flag)
                    { keepon = FALSE;
                      break;
                    }
               }
              if (keepon)           /* found it! mark u as nonprincipal */
               { rep[u] = v;
                 /* printf(" >> mapping %d onto %d\n", u, v); */
                 vtype[u] = 4;
                 u = next[u];
                 next[ulast] = u;  /* remove u from bin */
               }
              else                 /* failed */
               { ulast = u;
                 u = next[u];
               }
            }
           v = next[v];        /* no more variables can be absorbed by v */
           flag++;             /* clear tmp vector for next round */
         }
      }
   }
 
  /* --------------------
     free working storage
     -------------------- */
  free(tmp); free(bin);
  free(next); free(key);
}

◆ findIndNodes()

void findIndNodes	(	gelim_t *	Gelim,
		PORD_INT *	reachset,
		PORD_INT	nreach,
		PORD_INT *	bin,
		PORD_INT *	next,
		PORD_INT *	tmp,
		PORD_INT *	pflag )

Definition at line 639 of file gelim.c.

{ PORD_INT *xadj, *adjncy, *vwght, *len, *elen, *parent, *score;
  PORD_INT nvtx, chk, keepon, u, v, w, wlast, i, j, jstart, jstop, jstep, jj, jjstop;
  nvtx = Gelim->G->nvtx;
  xadj = Gelim->G->xadj;
  adjncy = Gelim->G->adjncy;
  vwght = Gelim->G->vwght;
  len = Gelim->len;
  elen = Gelim->elen;
  parent = Gelim->parent;
  score = Gelim->score;
 
#ifdef DEBUG
  printf("Checking reachset for indistinguishable variables\n");
#endif
 
  /* -----------------------------------------------------------------------
     compute checksums for all principal variables on reachset and fill bins
     NOTE: checksums are stored in parent vector
     ----------------------------------------------------------------------- */
  for (i = 0; i < nreach; i++)
   { u = reachset[i];
     chk = 0;
     jstart = xadj[u];
     jstop = jstart + len[u];
     /* Modified by JYL: 16 march 2005:
      * This code was failing in case of
      * overflow.
     for (j = jstart; j < jstop; j++)
         chk += adjncy[j];
     chk = chk % nvtx;
     */
     jstep=max(1000000000/nvtx,1);
     for (j = jstart; j < jstop; j+=jstep)
     {
       jjstop = min(jstop, j+jstep);
       for (jj = j; jj < jjstop; jj++)
         chk += adjncy[jj];
       chk = chk % nvtx;
     }
 
     parent[u] = chk;
     /* JYL: temporary:
        if (parent[u] < - 10)
        printf("Probleme %d \n",chk);*/
     next[u] = bin[chk];
     bin[chk] = u;
   }
 
  /* -----------------------
     supervariable detection
     ----------------------- */
  for (i = 0; i < nreach; i++)
   { u = reachset[i];
     if (vwght[u] > 0)         /* u is a principal variable */
      { chk = parent[u];       /* search bin[chk] for ind. nodes */
        v = bin[chk];          /* okay, v is the first node in this bin */
        bin[chk] = -1;         /* no further examinations of this bin */
        while (v != -1)
         { jstart = xadj[v];
           jstop = xadj[v] + len[v];
           for (j = jstart; j < jstop; j++)
             tmp[adjncy[j]] = *pflag;
           w = next[v];        /* v is principal and w is a potential */
           wlast = v;          /* nonprincipal variable */
           while (w != -1)
            { keepon = TRUE;
              if ((len[w] != len[v]) || (elen[w] != elen[v])
                 || ((score[w] < 0) && (score[v] >= 0))
                 || ((score[w] >= 0) && (score[v] < 0)))
                keepon = FALSE;
              if (keepon)
               { for (jj = xadj[w]; jj < xadj[w] + len[w]; jj++)
                   if (tmp[adjncy[jj]] < *pflag)
                    { keepon = FALSE;
                      break;
                    }
               }
              if (keepon)                /* found it! mark w as nonprincipal */
               { parent[w] = v;          /* representative of w is v */
                 /* Temporary JY
            if (parent[w] < - 10)
                   printf("Probleme\n");
                  */
#ifdef DEBUG
                 printf(" non-principal variable %d (score %d) mapped onto "
                        "%d (score %d)\n", w, score[w], v, score[v]); 
#endif
 
                 vwght[v] += vwght[w];   /* add weight of w */
                 vwght[w] = 0;
                 xadj[w] = -1;           /* w's adjacency list can be over- */
                 score[w] = -2;          /* written during next crunch */
                 w = next[w];
                 next[wlast] = w;        /* remove w from bin */
               }
              else                       /* failed */
               { wlast = w;
                 w = next[w];
               }
            }
           v = next[v];        /* no more variables can be absorbed by v */
           (*pflag)++;         /* clear tmp vector for next round */
         }
      }
   }
 
  /* -------------------------------------------------------
     re-initialize parent vector for all principal variables
     ------------------------------------------------------- */
  for (i = 0; i < nreach; i++)
   { u = reachset[i];
     if (vwght[u] > 0)
       parent[u] = -1;
   }
}   

◆ findPseudoPeripheralDomain()

PORD_INT findPseudoPeripheralDomain	(	domdec_t *	dd,
		PORD_INT	domain )

Definition at line 152 of file ddbisect.c.

{ PORD_INT *xadj, *adjncy, *vtype, *level, *queue;
  PORD_INT nvtx, qhead, qtail, nlev, lastdomain, u, v, i, istart, istop;
 
  nvtx = dd->G->nvtx;
  xadj = dd->G->xadj;
  adjncy = dd->G->adjncy;
  vtype = dd->vtype;
 
  /* ------------------------
     allocate working storage
     ------------------------ */
  mymalloc(level, nvtx, PORD_INT);
  mymalloc(queue, nvtx, PORD_INT);
 
  /* ---------------------------------------
     find a domain with maximal excentricity
     --------------------------------------- */
  nlev = 0; lastdomain = domain;
  while (TRUE)
   { for (u = 0; u < nvtx; u++)
       level[u] = -1;
     queue[0] = domain; level[domain] = 0;
     qhead = 0; qtail = 1;
     while (qhead != qtail)
      { u = queue[qhead++];
        if (vtype[u] == 1)     /* remember last domain */
          lastdomain = u;
        istart = xadj[u];
        istop = xadj[u+1];
        for (i = istart; i < istop; i++)
         { v = adjncy[i];
           if (level[v] == -1)
            { queue[qtail++] = v;
              level[v] = level[u] + 1;
            }
         }
      }
     if (level[lastdomain] > nlev)
      { nlev = level[lastdomain];
        domain = lastdomain;
      }
     else break;
   }
 
  /* -------------------------------
     free working storage and return
     ------------------------------- */
  free(level); free(queue);
  return(domain);
}

◆ firstPostorder()

PORD_INT firstPostorder ( elimtree_t * T )

Definition at line 213 of file tree.c.

{ PORD_INT *firstchild, J;
 
  firstchild = T->firstchild;
 
  if ((J = T->root) != -1)
    while (firstchild[J] != -1)
      J = firstchild[J];
  return(J);
}

◆ firstPostorder2()

PORD_INT firstPostorder2	(	elimtree_t *	T,
		PORD_INT	root )

Definition at line 228 of file tree.c.

{ PORD_INT *firstchild, J;
 
  firstchild = T->firstchild;
 
  if ((J = root) != -1)
    while (firstchild[J] != -1)
      J = firstchild[J];
  return(J);
}

◆ firstPreorder()

PORD_INT firstPreorder ( elimtree_t * T )

Definition at line 263 of file tree.c.

{
  return(T->root);
}

◆ forwardSubst1x1()

void forwardSubst1x1	(	factorMtx_t *	,
		FLOAT *	)

◆ forwardSubst1x1KernelPAR()

void forwardSubst1x1KernelPAR	(	topology_t *	,
		mask_t *	,
		PORD_INT	,
		PORD_INT	,
		factorMtx_t *	,
		FLOAT *	,
		FLOAT *	)

◆ forwardSubst1x1NEW()

void forwardSubst1x1NEW	(	factorMtx_t *	,
		FLOAT *	)

◆ forwardSubst1x1PAR()

void forwardSubst1x1PAR	(	topology_t *	,
		mask_t *	,
		mapping_t *	,
		factorMtx_t *	,
		FLOAT *	,
		FLOAT *	)

◆ freeBipartiteGraph()

void freeBipartiteGraph ( gbipart_t * Gbipart )

Definition at line 81 of file gbipart.c.

{
  freeGraph(Gbipart->G);
  free(Gbipart);
}

◆ freeBucket()

void freeBucket ( bucket_t * bucket )

Definition at line 78 of file bucket.c.

{
  free(bucket->bin);
  free(bucket->next);
  free(bucket->last);
  free(bucket->key);
  free(bucket);
}

◆ freeBuffer()

void freeBuffer ( buffer_t * )

◆ freeCSS()

void freeCSS ( css_t * css )

Definition at line 77 of file symbfac.c.

{
  free(css->xnzl);
  free(css->xnzlsub);
  if (css->owned)
    free(css->nzlsub);
  free(css);
}

◆ freeDenseMtx()

void freeDenseMtx ( denseMtx_t * )

◆ freeDomainDecomposition()

void freeDomainDecomposition ( domdec_t * dd )

Definition at line 120 of file ddcreate.c.

{
  freeGraph(dd->G);
  free(dd->vtype);
  free(dd->color);
  free(dd->map);
  free(dd);
}

◆ freeElimGraph()

void freeElimGraph ( gelim_t * Gelim )

Definition at line 128 of file gelim.c.

{
  freeGraph(Gelim->G);
  free(Gelim->len);
  free(Gelim->elen);
  free(Gelim->parent);
  free(Gelim->degree);
  free(Gelim->score);
  free(Gelim);
}

◆ freeElimTree()

void freeElimTree ( elimtree_t * T )

Definition at line 132 of file tree.c.

{
  free(T->ncolfactor);
  free(T->ncolupdate);
  free(T->parent);
  free(T->firstchild);
  free(T->silbings);
  free(T->vtx2front);
  free(T);
}

◆ freeFactorMtx()

void freeFactorMtx ( factorMtx_t * L )

Definition at line 465 of file symbfac.c.

{
  freeCSS(L->css);
  freeFrontSubscripts(L->frontsub);
  free(L->nzl);
  free(L->perm);
  free(L);
}

◆ freeFrontSubscripts()

void freeFrontSubscripts ( frontsub_t * frontsub )

Definition at line 286 of file symbfac.c.

{
  freeElimTree(frontsub->PTP);
  free(frontsub->xnzf);
  free(frontsub->nzfsub);
  free(frontsub);
}

◆ freeGbisect()

void freeGbisect ( gbisect_t * Gbisect )

Definition at line 77 of file gbisect.c.

{
  free(Gbisect->color);
  free(Gbisect);
}

◆ freeGraph()

void freeGraph ( graph_t * G )

Definition at line 73 of file graph.c.

{ 
  free(G->xadj);
  free(G->adjncy);
  free(G->vwght);
  free(G);
}

◆ freeInputMtx()

void freeInputMtx ( inputMtx_t * )

◆ freeMapping()

void freeMapping ( mapping_t * )

◆ freeMask()

void freeMask ( mask_t * )

◆ freeMinPriority()

void freeMinPriority ( minprior_t * minprior )

Definition at line 107 of file minpriority.c.

{
  freeElimGraph(minprior->Gelim);
  freeBucket(minprior->bucket);
  free(minprior->stageinfo);
  free(minprior->reachset);
  free(minprior->auxaux);
  free(minprior->auxbin);
  free(minprior->auxtmp);
  free(minprior);
}

◆ freeMultisector()

void freeMultisector ( multisector_t * ms )

Definition at line 86 of file multisector.c.

{
  free(ms->stage);
  free(ms);
}

◆ freeNDnode()

void freeNDnode ( nestdiss_t * nd )

Definition at line 96 of file nestdiss.c.

{
  free(nd->intvertex);
  free(nd->intcolor);
  free(nd);
}

◆ freeNDtree()

void freeNDtree ( nestdiss_t * ndroot )

Definition at line 260 of file nestdiss.c.

{ nestdiss_t *nd, *parent;
 
  /* ------------------------------------------------------
     to remove the nested dissection tree under root ndroot
     visit the nodes in post-order
     ------------------------------------------------------ */
  for (nd = ndroot; nd->childB != NULL; nd = nd->childB);
  while (nd != ndroot)
   { parent = nd->parent;
     if ((parent == NULL) || (parent->childB == NULL)
        || (parent->childW == NULL))
       { fprintf(stderr, "\nError in function removeNDtree\n"
              "  nested dissection tree corrupted\n");
         quit();
       }
     if (parent->childB == nd)  /* left subtree of parent visited */
      { freeNDnode(nd);         /* free root of left subtree and goto right */
        for (nd = parent->childW; nd->childB != NULL; nd = nd->childB);
      }
     else                       /* right subtree of parent visited */
      { freeNDnode(nd);         /* free root of right subtree and goto parent */
        nd = parent;
      }
   }
}

◆ freeTopology()

void freeTopology ( topology_t * )

◆ freeWorkspaceForDenseMtx()

void freeWorkspaceForDenseMtx ( workspace_t * )

◆ fundamentalFronts()

elimtree_t * fundamentalFronts ( elimtree_t * T )

Definition at line 553 of file tree.c.

{ elimtree_t *T2;
  PORD_INT        *ncolfactor, *ncolupdate, *parent, *firstchild, *silbings;
  PORD_INT        *frontmap, nfronts, cnfronts, J, child;
 
  nfronts = T->nfronts;
  ncolfactor = T->ncolfactor;
  ncolupdate = T->ncolupdate;
  parent = T->parent;
  firstchild = T->firstchild;
  silbings = T->silbings;
 
  /* -------------------------
     set up the working arrays
     ------------------------- */
  mymalloc(frontmap, nfronts, PORD_INT);
 
  /* -----------------------------
     search the fundamental fronts
     ----------------------------- */
  cnfronts = 0;
  J = T->root;
  while (J != -1)
   { while (firstchild[J] != -1)
       J = firstchild[J];
     frontmap[J] = cnfronts++;
     while ((silbings[J] == -1) && (parent[J] != -1))
      { J = parent[J];
        child = firstchild[J];
        if ((silbings[child] != -1)
            || (ncolupdate[child] != ncolfactor[J] + ncolupdate[J]))
          frontmap[J] = cnfronts++;
        else
          frontmap[J] = frontmap[child];
      }
     J = silbings[J];
   }
 
  /* ------------------------------
     construct new elimination tree
     ------------------------------ */
  T2 = compressElimTree(T, frontmap, cnfronts);
 
  /* ----------------------
     free memory and return
     ---------------------- */
  free(frontmap);
  return(T2);
}

◆ getWorkspaceForDenseMtx()

FLOAT * getWorkspaceForDenseMtx	(	workspace_t *	,
		PORD_INT	)

◆ greg_pord()

PORD_INT greg_pord	(	PORD_INT	,
		PORD_INT	,
		PORD_INT *	,
		PORD_INT *	,
		PORD_INT *	,
		PORD_INT *	,
		PORD_INT *	)

◆ improveDDSep()

void improveDDSep ( domdec_t * dd )

Definition at line 607 of file ddbisect.c.

{ bucket_t *b_bucket, *w_bucket;
  PORD_INT      *xadj, *adjncy, *vwght, *vtype, *color, *cwght;
  PORD_INT      *tmp_color, *deltaS, *deltaB, *deltaW;
  PORD_INT      nvtx, weight, tmp_S, tmp_B, tmp_W;
  PORD_INT      pos, bestglobalpos, badflips, b_domain, w_domain, domain, nxtdomain;
  PORD_INT      fhead, ftail, u, v, i, istart, istop;
  FLOAT    bestglobalvalue, b_value, w_value, value;
 
  /* ======================================================================
     vtype[u]: (u domain)
         1 => color of domain u has not been changed
       < 0 => points to next domain in flipping list
              (fhead points to first, ftail points to last domain in list)
       = 0 => domain is last domain in flipping list
     ====================================================================== */
 
  nvtx = dd->G->nvtx;
  xadj = dd->G->xadj;
  adjncy = dd->G->adjncy;
  vwght = dd->G->vwght;
  vtype = dd->vtype;
  color = dd->color;
  cwght = dd->cwght;
 
  mymalloc(tmp_color, nvtx, PORD_INT);
  mymalloc(deltaS, nvtx, PORD_INT);
  mymalloc(deltaB, nvtx, PORD_INT);
  mymalloc(deltaW, nvtx, PORD_INT);
 
OUTER_LOOP_START:
 
  /* ----------------------------------------------------------------------
     copy data of actual bisection and initialize buckets and flipping list
     ---------------------------------------------------------------------- */
  tmp_S = cwght[GRAY];
  tmp_B = cwght[BLACK];
  tmp_W = cwght[WHITE];
  bestglobalpos = badflips = 0;
  bestglobalvalue = F(tmp_S, tmp_B, tmp_W);
 
  b_bucket = setupBucket(nvtx, nvtx, (nvtx >> 1));
  w_bucket = setupBucket(nvtx, nvtx, (nvtx >> 1));
 
  fhead = 0; ftail = -1;
  pos = 0;
 
  /* ----------------------------------------------------------
     initialize tmp_color, deltaB, and deltaW for all multisecs
     ---------------------------------------------------------- */
  for (u = 0; u < nvtx; u++)
    if (vtype[u] == 2)
     { deltaB[u] = deltaW[u] = 0;
       istart = xadj[u];
       istop = xadj[u+1];
       for (i = istart; i < istop; i++)
        { v = adjncy[i];
          if (color[v] == BLACK) deltaB[u]++;
          else deltaW[u]++;
        }
       if ((deltaB[u] > 0) && (deltaW[u] > 0))  /* update multisec coloring */
         tmp_color[u] = GRAY;
       else if (deltaB[u] > 0) tmp_color[u] = BLACK;
       else tmp_color[u] = WHITE;
       color[u] = tmp_color[u];
     }
 
  /* -----------------------------------------------------------------
     initialize tmp_color, deltaS,B,W for all domains and fill buckets
     ----------------------------------------------------------------- */
  for (u = 0; u < nvtx; u++)
    if (vtype[u] == 1)
     { tmp_color[u] = color[u];
       if (tmp_color[u] == BLACK)          /* domain may be flipped to WHITE */
        { deltaW[u] = vwght[u]; deltaB[u] = -deltaW[u]; deltaS[u] = 0;
          istart = xadj[u];
          istop = xadj[u+1];
          for (i = istart; i < istop; i++)
           { v = adjncy[i];                /* tmp_color[v] e {GRAY, BLACK} */
             weight = vwght[v];
             if (tmp_color[v] == BLACK)    /* multisec v will move into S */
              { deltaB[u] -= weight;
                deltaS[u] += weight;
              }
             else if (deltaB[v] == 1)      /* multisec v will move into W */
              { deltaW[u] += weight;
                deltaS[u] -= weight;
                deltaB[v] = -(u+1);
              }
           }
          insertBucket(b_bucket, deltaS[u], u);
        }
       if (tmp_color[u] == WHITE)          /* domain may be flipped to BLACK */
        { deltaB[u] = vwght[u]; deltaW[u] = -deltaB[u]; deltaS[u] = 0;
          istart = xadj[u];
          istop = xadj[u+1];
          for (i = istart; i < istop; i++)
           { v = adjncy[i];                /* tmp_color[v] e {GRAY, WHITE} */
             weight = vwght[v];
             if (tmp_color[v] == WHITE)    /* multisec v will move into S */
              { deltaW[u] -= weight;
                deltaS[u] += weight;
              }
             else if (deltaW[v] == 1)      /* multisec v will move into B */
              { deltaB[u] += weight;
                deltaS[u] -= weight;
                deltaW[v] = -(u+1);
              }
           }
          insertBucket(w_bucket, deltaS[u], u);
        }
     }
 
#ifdef DEBUG
  printf("starting inner loop: b_bucket->nobj %d, w_bucket->nobj %d\n",
         b_bucket->nobj, w_bucket->nobj);
  waitkey();
#endif
 
INNER_LOOP_START:
 
  /* -------------------------------------------
     extract best domain from b_bucket, w_bucket
     ------------------------------------------- */
  b_value = w_value = MAX_FLOAT;
  if ((b_domain = minBucket(b_bucket)) != -1)
   { b_value = F((tmp_S+deltaS[b_domain]), (tmp_B+deltaB[b_domain]),
                 (tmp_W+deltaW[b_domain]));
 
#ifdef DEBUG
     printf("best black domain: %d, deltaS %d, deltaB %d, deltaW %d, "
            "cost %7.2f\n", b_domain, deltaS[b_domain], deltaB[b_domain],
            deltaW[b_domain], b_value);
#endif
   }
  if ((w_domain = minBucket(w_bucket)) != -1)
   { w_value = F((tmp_S+deltaS[w_domain]), (tmp_B+deltaB[w_domain]),
                 (tmp_W+deltaW[w_domain]));
 
#ifdef DEBUG
     printf("best white domain: %d, deltaS %d, deltaB %d, deltaW %d, "
            "cost %7.2f\n", w_domain, deltaS[w_domain], deltaB[w_domain],
            deltaW[w_domain], w_value);
#endif
   }
 
  if ((b_domain == ERR) && (w_domain == ERR)) goto INNER_LOOP_END;
 
  if (b_value + EPS < w_value)
   { domain = b_domain; value = b_value;
     removeBucket(b_bucket, domain);
   }
  else
   { domain = w_domain; value = w_value;
     removeBucket(w_bucket, domain);
   }
 
#ifdef DEBUG
  printf(" domain %d removed from bucket\n", domain);
#endif
 
  /* -------------------------------------------------------------------
     flip the color of domain and put it in list of log. flipped domains
     ------------------------------------------------------------------- */
  if (ftail != -1)
    vtype[ftail] = -(domain+1);   /* append domain */
  else fhead = -(domain+1);       /* list starts with domain */
  vtype[domain] = 0;              /* mark end of list */
  ftail = domain;                 /* domain is last element in list */
 
  if (tmp_color[domain] == BLACK)
   { tmp_color[domain] = WHITE;
     updateB2W(w_bucket,b_bucket,dd,domain,tmp_color,deltaW,deltaB,deltaS);
   }
  else if (tmp_color[domain] == WHITE)
   { tmp_color[domain] = BLACK;
     updateW2B(w_bucket,b_bucket,dd,domain,tmp_color,deltaW,deltaB,deltaS);
   }
  tmp_S += deltaS[domain];
  tmp_B += deltaB[domain];
  tmp_W += deltaW[domain];
 
  pos++;
  if (value + EPS < bestglobalvalue)
   { bestglobalvalue = value;
     bestglobalpos = pos;
     badflips = 0;
   }
  else badflips++;
  if (badflips < MAX_BAD_FLIPS) goto INNER_LOOP_START;
 
INNER_LOOP_END:
 
  /* --------------------------------------------
     end of inner loop: now do the physical flips
     -------------------------------------------- */
  pos = 0;
  nxtdomain = fhead;
  while (nxtdomain != 0)
   { domain = -nxtdomain - 1;
     if (pos < bestglobalpos)
      { if (color[domain] == BLACK) color[domain] = WHITE;
        else color[domain] = BLACK;
        cwght[GRAY] += deltaS[domain];
        cwght[BLACK] += deltaB[domain];
        cwght[WHITE] += deltaW[domain];
        pos++;
      }
     nxtdomain = vtype[domain];
     vtype[domain] = 1;
   }
 
  /* ----------------------------------------------
     partition improved => re-start the whole stuff
     ---------------------------------------------- */
#ifdef DEBUG
  printf(" INNER_LOOP_END (#pyhs. flips %d): S %d, B %d, W %d (%7.2f)\n",
         bestglobalpos, cwght[GRAY], cwght[BLACK], cwght[WHITE],
         bestglobalvalue);
  waitkey();
#endif
 
  /* JY: moved next instruction after the two
   *     freeBucket instructions because
   *     this was the cause of a memory leak.
   * if (bestglobalpos > 0) goto OUTER_LOOP_START;
   */
 
  freeBucket(b_bucket);
  freeBucket(w_bucket);
 
  if (bestglobalpos > 0) goto OUTER_LOOP_START;
  free(tmp_color); free(deltaS); free(deltaB); free(deltaW);
}

◆ initFactorMtx()

void initFactorMtx	(	factorMtx_t *	L,
		inputMtx_t *	PAP )

Definition at line 510 of file symbfac.c.

{ elimtree_t *PTP;
  frontsub_t *frontsub;
  css_t      *css;
  PORD_INT        *ncolfactor;
  FLOAT      *nzl, *nza, *diag;
  PORD_INT        *xnzl, *nzlsub, *xnzlsub, *xnza, *nzasub, *xnzf, *nzfsub;
  PORD_INT        nelem, K, k, kstart, h, hstart, dis, i, istart, istop;
  PORD_INT        firstcol, lastcol;
 
  nelem = L->nelem; 
  nzl = L->nzl;
  css = L->css;
  xnzl = css->xnzl;
  nzlsub = css->nzlsub;
  xnzlsub = css->xnzlsub;
 
  frontsub = L->frontsub;
  PTP = frontsub->PTP;
  ncolfactor = PTP->ncolfactor;
  xnzf = frontsub->xnzf;
  nzfsub = frontsub->nzfsub;
 
  diag = PAP->diag;
  nza = PAP->nza;
  xnza = PAP->xnza;
  nzasub = PAP->nzasub;
 
  /* ------------------------------------
     set all numerical values of L to 0.0
     ------------------------------------ */
  for (i = 0; i < nelem; i++)
    nzl[i] = 0.0;
 
  /* --------------------------------------------
     init. factor matrix with the nonzeros of PAP
     -------------------------------------------- */
  for (K = firstPostorder(PTP); K != -1; K = nextPostorder(PTP, K))
   { firstcol = nzfsub[xnzf[K]];
     lastcol = firstcol + ncolfactor[K];
     for (k = firstcol; k < lastcol; k++)
      { istart = xnza[k];
        istop = xnza[k+1];
        kstart = xnzl[k];
        hstart = xnzlsub[k];
        h = hstart;
        for (i = istart; i < istop; i++)
         { for (; nzlsub[h] != nzasub[i]; h++);
           dis = h - hstart;
           nzl[kstart+dis] = nza[i];
         }
        nzl[kstart] = diag[k];
      }
   }
}

◆ initFactorMtxNEW()

void initFactorMtxNEW	(	factorMtx_t *	L,
		inputMtx_t *	PAP )

Definition at line 570 of file symbfac.c.

{ elimtree_t *PTP;
  frontsub_t *frontsub;
  css_t      *css;
  PORD_INT        *ncolfactor;
  FLOAT      *nzl, *nza, *diag, *entriesL;
  PORD_INT        *xnzl, *xnza, *nzasub, *xnzf, *nzfsub;
  PORD_INT        *tmp, neqs, nelem, K, k, len, row, i, istart, istop;
  PORD_INT        firstcol, lastcol;
 
  nelem = L->nelem;
  nzl = L->nzl;
  css = L->css;
  xnzl = css->xnzl;
 
  frontsub = L->frontsub;
  PTP = frontsub->PTP;
  ncolfactor = PTP->ncolfactor;
  xnzf = frontsub->xnzf;
  nzfsub = frontsub->nzfsub;
 
  neqs = PAP->neqs;
  diag = PAP->diag;
  nza = PAP->nza;
  xnza = PAP->xnza;
  nzasub = PAP->nzasub;
 
  /* ------------------------
     allocate working storage
     ------------------------ */
  mymalloc(tmp, neqs, PORD_INT);
 
  /* ------------------------------------
     set all numerical values of L to 0.0
     ------------------------------------ */
  for (i = 0; i < nelem; i++)
    nzl[i] = 0.0;
 
  /* --------------------------------------------
     init. factor matrix with the nonzeros of PAP
     -------------------------------------------- */
  for (K = firstPostorder(PTP); K != -1; K = nextPostorder(PTP, K))
   { len = 0;
     istart = xnzf[K];
     istop = xnzf[K+1];
     for (i = istart; i < istop; i++)
       tmp[nzfsub[i]] = len++;
 
     firstcol = nzfsub[istart];
     lastcol = firstcol + ncolfactor[K];
     entriesL = nzl + xnzl[firstcol];
     for (k = firstcol; k < lastcol; k++)
      { istart = xnza[k];
        istop = xnza[k+1];
        for (i = istart; i < istop; i++)
         { row = nzasub[i];
           entriesL[tmp[row]] = nza[i];
         }
        entriesL[tmp[k]] = diag[k];
        entriesL += --len;
      }
   }
 
  /* --------------------
     free working storage
     -------------------- */
  free(tmp);
}

◆ initFactorMtxPAR()

void initFactorMtxPAR	(	mask_t *	,
		mapping_t *	,
		factorMtx_t *	,
		inputMtx_t *	)

◆ initFchSilbRoot()

void initFchSilbRoot ( elimtree_t * T )

Definition at line 414 of file tree.c.

{ PORD_INT *parent, *firstchild, *silbings, nfronts, J, pJ;
 
  nfronts = T->nfronts;
  parent = T->parent;
  firstchild = T->firstchild;
  silbings = T->silbings;
 
  for (J = 0; J < nfronts; J++)
    silbings[J] = firstchild[J] = -1;
 
  for (J = nfronts-1; J >= 0; J--)
    if ((pJ = parent[J]) != -1)
     { silbings[J] = firstchild[pJ];
       firstchild[pJ] = J;
     }
    else
     { silbings[J] = T->root;
       T->root = J;
     }
}

◆ initialDDSep()

void initialDDSep ( domdec_t * dd )

Definition at line 348 of file ddbisect.c.

{  PORD_INT *vtype, *color, *cwght;
   PORD_INT nvtx, totvwght, domain, u;
 
   nvtx = dd->G->nvtx;
   totvwght = dd->G->totvwght;
   vtype = dd->vtype;
   color = dd->color;
   cwght = dd->cwght;
 
  /* --------------------------------------------------------
     initializations (all nodes are colored white by default)
     -------------------------------------------------------- */
  cwght[GRAY] = 0;
  cwght[BLACK] = 0;
  cwght[WHITE] = totvwght;
  for (u = 0; u < nvtx; u++)
    color[u] = WHITE;
 
  /* ----------------------------------------------------------------------
     scan over connected components and create level based vertex separator
     ---------------------------------------------------------------------- */
  for (u = 0; u < nvtx; u++)
    if ((vtype[u] == 1) && (color[u] == WHITE))
     { domain = findPseudoPeripheralDomain(dd, u);
       constructLevelSep(dd, domain);
       if (cwght[BLACK] >= cwght[WHITE])
         break;
     }
}

◆ initialDomainDecomposition()

domdec_t * initialDomainDecomposition	(	graph_t *	G,
		PORD_INT *	map,
		PORD_INT *	vtype,
		PORD_INT *	rep )

Definition at line 366 of file ddcreate.c.

{ domdec_t *dd;
  PORD_INT      *xadj, *adjncy, *vwght, *xadjdd, *adjncydd, *vwghtdd, *vtypedd;
  PORD_INT      *tmp, *bin, nvtx, nedges, nvtxdd, nedgesdd, ndom, domwght, flag;
  PORD_INT      i, j, jstart, jstop, u, v, w;
 
  nvtx = G->nvtx;
  nedges = G->nedges;
  xadj = G->xadj;
  adjncy = G->adjncy;
  vwght = G->vwght;
 
  /* ------------------------
     allocate working storage
     ------------------------ */
  mymalloc(tmp, nvtx, PORD_INT);
  mymalloc(bin, nvtx, PORD_INT);
  for (u = 0; u < nvtx; u++)
   { tmp[u] = -1;
     bin[u] = -1;
   }
 
  /* -------------------------------------------------------------
     allocate memory for the dd using upper bounds nvtx and nedges
     ------------------------------------------------------------- */
  dd = newDomainDecomposition(nvtx, nedges);
  xadjdd = dd->G->xadj;
  adjncydd = dd->G->adjncy;
  vwghtdd = dd->G->vwght;
  vtypedd = dd->vtype;
 
  /* -------------------------------------------------------
     put all nodes u belonging to representative v in bin[v]
     ------------------------------------------------------- */
  for (u = 0; u < nvtx; u++)
   { v = rep[u];
     if (u != v)
      { bin[u] = bin[v];
        bin[v] = u;
      }
   }
 
  /* ----------------------------------------------
     and now build the initial domain decomposition
     ---------------------------------------------- */
  flag = 1;
  nedgesdd = nvtxdd = 0;
  ndom = domwght = 0;
  for (u = 0; u < nvtx; u++)
    if (rep[u] == u)
     { xadjdd[nvtxdd] = nedgesdd;
       vtypedd[nvtxdd] = vtype[u];
       vwghtdd[nvtxdd] = 0;
       tmp[u] = flag;
 
       /* find all cluster that are adjacent to u in dom. dec. */
       v = u;
       do
        { map[v] = nvtxdd;
          vwghtdd[nvtxdd] += vwght[v];
          jstart = xadj[v];
          jstop = xadj[v+1];
          for (j = jstart; j < jstop; j++)
           { w = adjncy[j];
             if ((vtype[w] != vtype[u]) && (tmp[rep[w]] != flag))
              { tmp[rep[w]] = flag;
                adjncydd[nedgesdd++] = rep[w];
              }
           }
          v = bin[v];
        } while (v != -1);
 
       if (vtypedd[nvtxdd] == 1)
        { ndom++;
          domwght += vwghtdd[nvtxdd];
        }
       nvtxdd++;
       flag++;
     }
 
  /* --------------------------------------------
     finalize the new domain decomposition object
     -------------------------------------------- */
  xadjdd[nvtxdd] = nedgesdd;
  dd->G->nvtx = nvtxdd;
  dd->G->nedges = nedgesdd;
  dd->G->type = WEIGHTED;
  dd->G->totvwght = G->totvwght;
  for (i = 0; i < nedgesdd; i++)
    adjncydd[i] = map[adjncydd[i]];
  for (u = 0; u < nvtxdd; u++)
    dd->color[u] = dd->map[u] = -1;
  dd->ndom = ndom;
  dd->domwght = domwght;
 
  /* -------------------------------
     free working storage and return
     ------------------------------- */
  free(tmp); free(bin);
  return(dd);
}

◆ initLocalIndices()

void initLocalIndices	(	denseMtx_t *	,
		PORD_INT *	,
		PORD_INT *	)

◆ initLocalIndicesPAR()

void initLocalIndicesPAR	(	denseMtx_t *	,
		PORD_INT *	,
		PORD_INT *	)

◆ initWorkspaceForDenseMtx()

workspace_t * initWorkspaceForDenseMtx	(	PORD_INT	,
		PORD_INT	)

◆ insertBucket()

void insertBucket	(	bucket_t *	bucket,
		PORD_INT	k,
		PORD_INT	item )

Definition at line 163 of file bucket.c.

{ PORD_INT s, nextitem;
 
  /* ------------------------------------
     check whether there are any problems
     ------------------------------------ */
  if (abs(k) >= MAX_INT - bucket->offset - 1)
   { fprintf(stderr, "\nError in function insertBucket\n"
          "  key %d too large/small for bucket\n", k);
     quit();
   }
  if (item > bucket->maxitem)
   { fprintf(stderr, "\nError in function insertBucket\n"
          "  item %d too large for bucket (maxitem is %d)\n", item,
          bucket->maxitem);
     quit();
   }
  if (bucket->key[item] != MAX_INT)
   { fprintf(stderr, "\nError in function insertBucket\n"
          "  item %d already in bucket\n", item);
     quit();
   }
 
  /* -------------------------------------
     determine the bin that holds the item
     ------------------------------------- */
  s = max(0, (k + bucket->offset));
  s = min(s, bucket->maxbin);
 
  /* --------------------------------------------------------------
     adjust minbin, increase nobj, and mark item as being in bucket
     -------------------------------------------------------------- */
  bucket->minbin = min(bucket->minbin, s);
  bucket->nobj++;
  bucket->key[item] = k;
 
  /* -----------------------------
     finally, insert item in bin s
     ----------------------------- */
  nextitem = bucket->bin[s];
  if (nextitem != -1)
    bucket->last[nextitem] = item;
  bucket->next[item] = nextitem;
  bucket->last[item] = -1;
  bucket->bin[s] = item;
}

◆ insertDownIntsWithStaticFloatKeys()

void insertDownIntsWithStaticFloatKeys	(	PORD_INT	n,
		PORD_INT *	array,
		FLOAT *	key )

Definition at line 58 of file sort.c.

{ PORD_INT   i, j, e;
  FLOAT ke;
 
  for (i = 1; i < n; i++)
   { e = array[i]; ke = key[e]; j = i;
     while ((j > 0) && (key[array[j-1]] < ke))
      { array[j] = array[j-1];
        j--;
      }
     array[j] = e;
   }
}

◆ insertUpFloatsWithIntKeys()

void insertUpFloatsWithIntKeys	(	PORD_INT	n,
		FLOAT *	array,
		PORD_INT *	key )

Definition at line 77 of file sort.c.

{ PORD_INT   i, j, ke;
  FLOAT e;
 
  for (i = 1; i < n; i++)
   { e = array[i]; ke = key[i]; j = i;
     while ((j > 0) && (key[j-1] > ke))
      { array[j] = array[j-1];
        key[j] = key[j-1];
        j--;
      }
     array[j] = e;
     key[j] = ke;
   }
}

◆ insertUpInts()

void insertUpInts	(	PORD_INT	n,
		PORD_INT *	array )

Definition at line 21 of file sort.c.

{ PORD_INT i, j, v;
 
  for (i = 1; i < n; i++)
   { v = array[i]; j = i;
     while ((j > 0) && (array[j-1] > v))
      { array[j] = array[j-1];
        j--;
      }
     array[j] = v;
   }
}

◆ insertUpIntsWithStaticIntKeys()

void insertUpIntsWithStaticIntKeys	(	PORD_INT	n,
		PORD_INT *	array,
		PORD_INT *	key )

Definition at line 39 of file sort.c.

{ PORD_INT   i, j, ke;
  PORD_INT   e;
 
  for (i = 1; i < n; i++)
   { e = array[i]; ke = key[e]; j = i;
     while ((j > 0) && (key[array[j-1]] > ke))
      { array[j] = array[j-1];
        j--;
      }
     array[j] = e;
   }
}

◆ justifyFronts()

PORD_INT justifyFronts ( elimtree_t * T )

Definition at line 741 of file tree.c.

{ PORD_INT *ncolfactor, *ncolupdate, *firstchild, *silbings, *minWspace, *list;
  PORD_INT nfronts, K, ncolfrontK, frontsizeK, wspace, child, nxtchild;
  PORD_INT count, m, s, i;
 
  nfronts = T->nfronts;
  ncolfactor = T->ncolfactor;
  ncolupdate = T->ncolupdate;
  firstchild = T->firstchild;
  silbings = T->silbings;
 
  /* -------------------------
     set up the working arrays
     ------------------------- */
  mymalloc(minWspace, nfronts, PORD_INT);
  mymalloc(list, nfronts, PORD_INT);
 
  /* ---------------------------------------------------------
     postorder traversal of the elimination tree to obtain the
     optimal justification of the children of each front
     ---------------------------------------------------------- */
  wspace = 0;
  for (K = firstPostorder(T); K != -1; K = nextPostorder(T, K))
   { ncolfrontK = ncolfactor[K] + ncolupdate[K];
     frontsizeK = (ncolfrontK * (ncolfrontK + 1)) >> 1;
 
     if ((child = firstchild[K]) == -1)
       minWspace[K] = frontsizeK;
     else
      { count = 0;
 
        /* sort children according to their minWspace value */
        while (child != -1)
         { list[count++] = child;
           child = silbings[child];
         }
        insertUpIntsWithStaticIntKeys(count, list, minWspace);
        firstchild[K] = -1;
        for (i = 0; i < count; i++)
         { child = list[i];
           silbings[child] = firstchild[K];
           firstchild[K] = child;
         }
 
        /* compute minWspace[K] */
        child = firstchild[K];
        nxtchild = silbings[child];
        m = s = minWspace[child];
        while (nxtchild != -1)
         { s = s - minWspace[child]
               + ((ncolupdate[child] * (ncolupdate[child] + 1)) >> 1)
               + minWspace[nxtchild];
           m = max(m, s);
           child = nxtchild;
           nxtchild = silbings[nxtchild];
         }
        s = s - minWspace[child]
            + ((ncolupdate[child] * (ncolupdate[child] + 1)) >> 1)
            + frontsizeK;
        minWspace[K] = max(m, s);
      }
 
     wspace = max(wspace, minWspace[K]);
   }
 
  /* ----------------------
     free memory and return
     ---------------------- */
  free(minWspace); free(list);
  return(wspace);
}

◆ listing()

void listing	(	mapping_t *	,
		PORD_INT	,
		PORD_INT	,
		PORD_INT	,
		FLOAT *	,
		FLOAT *	)

◆ maximumFlow()

void maximumFlow	(	gbipart_t *	Gbipart,
		PORD_INT *	flow,
		PORD_INT *	rc )

Definition at line 322 of file gbipart.c.

{ PORD_INT *xadj, *adjncy, *vwght, *parent, *marker, *queue;
  PORD_INT nedges, u, v, x, y, nX, nY, j, i, istart, istop;
  PORD_INT qhead, qtail, capacity;
 
  nedges = Gbipart->G->nedges;
  xadj = Gbipart->G->xadj;
  adjncy = Gbipart->G->adjncy;
  vwght = Gbipart->G->vwght;
  nX = Gbipart->nX;
  nY = Gbipart->nY;
 
  mymalloc(parent, (nX+nY), PORD_INT);
  mymalloc(marker, (nX+nY), PORD_INT);
  mymalloc(queue, (nX+nY), PORD_INT);
 
  /* -------------------------------------
     initialize flow and residual capacity
     ------------------------------------- */
  for (u = 0; u < nX+nY; u++)
    rc[u] = vwght[u];
  for (i = 0; i < nedges; i++)
    flow[i] = 0;
 
  /* --------------------------------------------------
     determine an initial flow in the bipartite network
     -------------------------------------------------- */
  for (x = 0; x < nX; x++)
   { istart = xadj[x];
     istop = xadj[x+1];
     for (i = istart; i < istop; i++)
      { y = adjncy[i];
        capacity = min(rc[x], rc[y]);
        if (capacity > 0)
         { rc[x] -= capacity;
           rc[y] -= capacity;
           flow[i] = capacity;
           for (j = xadj[y]; adjncy[j] != x; j++);
           flow[j] = -capacity;
         }
        if (rc[x] == 0) break;
      }
   }
 
  /* -----------------------------------------------------------
     construct maximum flow in bipartite network (Edmonds, Karp)
     ----------------------------------------------------------- */
  while (TRUE)
   { for (u = 0; u < nX+nY; u++)
       parent[u] = marker[u] = -1;
     qhead = qtail = 0;                /* fill queue with free X nodes */
     for (x = 0; x < nX; x++)
       if (rc[x] > 0)
        { queue[qtail++] = x;
          parent[x] = x;
        }
 
      /* ---------------------------------------------------------
         breadth first search to find the shortest augmenting path
         --------------------------------------------------------- */
      capacity = 0;
      while (qhead != qtail)
       { u = queue[qhead++];
         istart = xadj[u];
         istop = xadj[u+1];
         for (i = istart; i < istop; i++)
          { v = adjncy[i];
            if ((parent[v] == -1) && ((v >= nX) || (flow[i] < 0)))
            /* v >= nX => u->v is a forward edge having infty capacity     */
            /* otherwise  u<-v is a backward edge and (v,u) must have      */
            /*            positive capacity (i.e. (u,v) has neg. capacity) */
             { parent[v] = u;
               marker[v] = i;
               queue[qtail++] = v;
               if ((v >= nX) && (rc[v] > 0))  /* found it! */
                { u = v;                      /* (v,sink) is below capacity */
                  capacity = rc[u];
                  while (parent[u] != u)      /* get minimal residual capa. */
                   { i = marker[u];
                     u = parent[u];
                     if (u >= nX)
                       capacity = min(capacity, -flow[i]);
                   }
                  capacity = min(capacity, rc[u]);
                  rc[v] -= capacity;          /* augment flow by min. rc */
                  while (parent[v] != v)
                   { i = marker[v];
                     u = parent[v];
                     flow[i] += capacity;
                     for (j = xadj[v]; adjncy[j] != u; j++);
                     flow[j] = -flow[i];
                     v = u;
                   }
                  rc[v] -= capacity;
                  qhead = qtail;              /* escape inner while loop */
                  break;
                }
             }
          }
       }
 
     if (capacity == 0)
       break;
   }
 
  free(parent); free(marker);
  free(queue);
}

◆ maximumMatching()

void maximumMatching	(	gbipart_t *	Gbipart,
		PORD_INT *	matching )

Definition at line 195 of file gbipart.c.

{ PORD_INT *xadj, *adjncy, *level, *marker, *queue, *stack;
  PORD_INT top, top2, u, x, x2, y, y2, nX, nY, i, istart, istop;
  PORD_INT qhead, qtail, max_level;
 
  xadj = Gbipart->G->xadj;
  adjncy = Gbipart->G->adjncy;
  nX = Gbipart->nX;
  nY = Gbipart->nY;
 
  mymalloc(level, (nX+nY), PORD_INT);
  mymalloc(marker, (nX+nY), PORD_INT);
  mymalloc(queue, nX, PORD_INT);
  mymalloc(stack, nY, PORD_INT);
 
  /* -------------------
     initialize matching
     ------------------- */
  for (u = 0; u < nX+nY; u++)
    matching[u] = FREE;
 
  /* ---------------------------------------------------
     construct maximal matching in bipartite graph (X,Y)
     --------------------------------------------------- */
  for (x = 0; x < nX; x++)
   { istart = xadj[x];
     istop = xadj[x+1];
     for (i = istart; i < istop; i++)
      { y = adjncy[i];
        if (matching[y] == FREE)
         { matching[x] = y;
           matching[y] = x;
           break;
         }
      }
   }
 
  /* --------------------------------------------------------------------
     construct maximum matching in bipartite graph (X,Y) (Hopcroft, Karp)
     -------------------------------------------------------------------- */
  while (TRUE)
   { for (u = 0; u < nX+nY; u++)
       level[u] = marker[u] = -1;
     qhead = qtail = 0;                /* fill queue with free X nodes */
     for (x = 0; x < nX; x++)
       if (matching[x] == FREE)
        { queue[qtail++] = x;
          level[x] = 0;
        }
 
     /* --------------------------------------------------------------
        breadth first search to construct layer network containing all
        vertex disjoint augmenting paths of minimal length
        -------------------------------------------------------------- */
     top = 0;
     max_level = MAX_INT;
     while (qhead != qtail)
      { x = queue[qhead++];                  /* note: queue contains only */
        if (level[x] < max_level)            /*       nodes from X        */
         { istart = xadj[x];
           istop = xadj[x+1];
           for (i = istart; i < istop; i++)
            { y = adjncy[i];
              if (level[y] == -1)
               { level[y] = level[x] + 1;
                 if (matching[y] == FREE)
                  { max_level = level[y];    /* note: stack contains only */
                    stack[top++] = y;        /*       nodes form Y        */
                  }
                 else if (level[y] < max_level)
                  { x2 = matching[y];
                    level[x2] = level[y] + 1;
                    queue[qtail++] = x2;
                  }
               }
            }
         }
      }
     if (top == 0) break;              /* no augmenting path found */ 
 
     /* ------------------------------------------------------------
        restricted depth first search to construct maximal number of
        vertex disjoint augmenting paths in layer network
        ------------------------------------------------------------ */
     while (top > 0)
      { top2 = top--;
        y = stack[top2-1];             /* get the next exposed node in Y */
        marker[y] = xadj[y];           /* points to next neighbor of y */
 
        while (top2 > top)
         { y = stack[top2-1];
           i = marker[y]++;
           if (i < xadj[y+1])          /* not all neighbors of y visited */
            { x = adjncy[i];
              if ((marker[x] == -1) && (level[x] == level[y]-1))
               { marker[x] = 0;
                 if (level[x] == 0)        /* augmenting path found */
                   while (top2 > top)      /* pop stack */
                    { y2 = stack[--top2];
                      x2 = matching[y2];   /*    / o == o        */
                      matching[x] = y2;    /*   /                */
                      matching[y2] = x;    /* x -- y2 == x2 -- y */
                      x = x2;              /*   \                */
                    }                      /*    \ o == o        */
                 else
                  { y2 = matching[x];
                    stack[top2++] = y2;
                    marker[y2] = xadj[y2];
                  }
               }
            }
           else top2--;
         }
      }
   }
 
  /* -------------------------------
     free working storage and return
     ------------------------------- */
  free(level); free(marker);
  free(queue); free(stack);
}

◆ mergeFronts()

elimtree_t * mergeFronts	(	elimtree_t *	T,
		PORD_INT	maxzeros )

Definition at line 607 of file tree.c.

{ elimtree_t *T2;
  PORD_INT        *ncolfactor, *ncolupdate, *firstchild, *silbings;
  PORD_INT        *frontmap, *newncolfactor, *nzeros, *rep;
  PORD_INT        nfronts, cnfronts, K, ncolfrontK, J, Jall, cost;
 
  nfronts = T->nfronts;
  ncolfactor = T->ncolfactor;
  ncolupdate = T->ncolupdate;
  firstchild = T->firstchild;
  silbings = T->silbings;
 
  /* -------------------------
     set up the working arrays
     ------------------------- */
  mymalloc(frontmap, nfronts, PORD_INT);
  mymalloc(newncolfactor, nfronts, PORD_INT);
  mymalloc(nzeros, nfronts, PORD_INT);
  mymalloc(rep, nfronts, PORD_INT);
  for (K = 0; K < nfronts; K++)
   { newncolfactor[K] = ncolfactor[K];
     nzeros[K] = 0;
     rep[K] = K;
   }
 
  /* -----------------------------------------------------
     perform a postorder traversal of the elimination tree
     ----------------------------------------------------- */
  for (K = firstPostorder(T); K != -1; K = nextPostorder(T, K))
    if (firstchild[K] != -1)
     { ncolfrontK = newncolfactor[K] + ncolupdate[K];
       Jall = 0;
       cost = 0;
       for (J = firstchild[K]; J != -1; J = silbings[J])
        { Jall += newncolfactor[J];
          cost -= newncolfactor[J] * newncolfactor[J];
          cost += 2*newncolfactor[J] * (ncolfrontK - ncolupdate[J]);
          cost += 2*nzeros[J];
        }
       cost += Jall * Jall;
       cost = cost / 2;
       if (cost < maxzeros)
        { for (J = firstchild[K]; J != -1; J = silbings[J])
           { rep[J] = K;
             newncolfactor[K] += newncolfactor[J];
           }
          nzeros[K] = cost;
        }
     }
 
  /* ----------------------------------
     construct frontmap from vector rep
     ---------------------------------- */
  cnfronts = 0;
  for (K = 0; K < nfronts; K++)
    if (rep[K] == K)
      frontmap[K] = cnfronts++;
    else
     { for (J = K; rep[J] != J; J = rep[J]);
       rep[K] = J;
     }
  for (K = 0; K < nfronts; K++)
    if ((J = rep[K]) != K)
      frontmap[K] = frontmap[J];
 
  /* ------------------------------
     construct new elimination tree
     ------------------------------ */
  T2 = compressElimTree(T, frontmap, cnfronts);
 
  /* ----------------------
     free memory and return
     ---------------------- */
  free(frontmap); free(newncolfactor);
  free(nzeros); free(rep);
  return(T2);
}

◆ mergeMultisecs()

void mergeMultisecs	(	graph_t *	G,
		PORD_INT *	vtype,
		PORD_INT *	rep )

Definition at line 280 of file ddcreate.c.

{ PORD_INT *xadj, *adjncy, *tmp, *queue;
  PORD_INT nvtx, qhead, qtail, flag, keepon, u, v, w, x;
  PORD_INT i, istart, istop, j, jstart, jstop;
 
  nvtx = G->nvtx;
  xadj = G->xadj;
  adjncy = G->adjncy;
 
  /* ------------------------
     allocate working storage
     ------------------------ */
  mymalloc(tmp, nvtx, PORD_INT);
  mymalloc(queue, nvtx, PORD_INT);
  for (u = 0; u < nvtx; u++)
    tmp[u] = -1;
 
  /* -------------------------------------------------------
     merge all adjacent multisecs that do not share a domain
     ------------------------------------------------------- */
  flag = 1;
  for (u = 0; u < nvtx; u++)
    if (vtype[u] == 2)
     { qhead = 0; qtail = 1;
       queue[0] = u;
       vtype[u] = -2;
 
       /* multisec u is the seed of a new cluster, mark all adj. domains */
       istart = xadj[u];
       istop = xadj[u+1];
       for (i = istart; i < istop; i++)
        { v = adjncy[i];
          if (vtype[v] == 1)
            tmp[rep[v]] = flag;
        }
 
       /* and now build the cluster */
       while (qhead != qtail)
        { v = queue[qhead++];
          istart = xadj[v];
          istop = xadj[v+1];
          for (i = istart; i < istop; i++)
           { keepon = TRUE;
             w = adjncy[i];
             if (vtype[w] == 2)
              { jstart = xadj[w];
                jstop = xadj[w+1];
                for (j = jstart; j < jstop; j++)
                 { x = adjncy[j];
                   if ((vtype[x] == 1) && (tmp[rep[x]] == flag))
                    { keepon = FALSE;
                      break;
                    }
                 }
                if (keepon)
                /* multisecs v and w have no domain in common; mark */
                /* all domains adjacent to w and put w in cluster u */
                 { for (j = jstart; j < jstop; j++)
                    { x = adjncy[j];
                      if (vtype[x] == 1) tmp[rep[x]] = flag;
                     }
                    queue[qtail++] = w;
                    rep[w] = u;
                    vtype[w] = -2;
                 }
              }
           }
        }
 
       /* clear tmp vector for next round */
       flag++;
     }
 
  /* ------------------------------------
     reset vtype and free working storage
     ------------------------------------ */
  for (u = 0; u < nvtx; u++)
    if (vtype[u] == -2)
      vtype[u] = 2;
  free(tmp); free(queue);
}

◆ minBucket()

PORD_INT minBucket ( bucket_t * bucket )

Definition at line 117 of file bucket.c.

{ PORD_INT *bin, *next, *key, maxbin, minbin, nobj;
  PORD_INT item, bestitem, bestkey;
 
  maxbin = bucket->maxbin;
  nobj = bucket->nobj;
  minbin = bucket->minbin;
  bin = bucket->bin;
  next = bucket->next;
  key = bucket->key;
 
  if (nobj > 0)
   { /* ---------------------------------------------
        get the first item from leftmost nonempty bin
        --------------------------------------------- */
     while (bin[minbin] == -1) minbin++;
     bucket->minbin = minbin;
     bestitem = bin[minbin];
     bestkey = minbin;
 
     /* --------------------------------------------------
        items in bins 0 and maxbin can have different keys
        => search for item with smallest key
        -------------------------------------------------- */
     if ((minbin == 0) || (minbin == maxbin))
      { item = next[bestitem];
        while (item != -1)
         { if (key[item] < bestkey)
            { bestitem = item;
              bestkey = key[item];
            }
           item = next[item];
         }
      }
     /* ---------------------------------
        return the item with smallest key
        --------------------------------- */
     return(bestitem);
   }
  else return(-1);
}

◆ myrank()

PORD_INT myrank ( void )

◆ newBipartiteGraph()

gbipart_t * newBipartiteGraph	(	PORD_INT	nX,
		PORD_INT	nY,
		PORD_INT	nedges )

Definition at line 66 of file gbipart.c.

{ gbipart_t *Gbipart;
 
  mymalloc(Gbipart, 1, gbipart_t);
  Gbipart->G = newGraph(nX+nY, nedges);
  Gbipart->nX = nX;
  Gbipart->nY = nY;
 
  return(Gbipart);
}

◆ newBucket()

bucket_t * newBucket	(	PORD_INT	maxbin,
		PORD_INT	maxitem,
		PORD_INT	offset )

Definition at line 56 of file bucket.c.

{ bucket_t *bucket;
 
  mymalloc(bucket, 1, bucket_t);
  mymalloc(bucket->bin, (maxbin+1), PORD_INT);
  mymalloc(bucket->next, (maxitem+1), PORD_INT);
  mymalloc(bucket->last, (maxitem+1), PORD_INT);
  mymalloc(bucket->key, (maxitem+1), PORD_INT);
 
  bucket->maxbin = maxbin;
  bucket->maxitem = maxitem;
  bucket->offset = offset;
  bucket->nobj = 0;
  bucket->minbin = MAX_INT;
 
  return(bucket);
}

◆ newBuffer()

buffer_t * newBuffer ( size_t )

◆ newCSS()

css_t * newCSS	(	PORD_INT	neqs,
		PORD_INT	nind,
		PORD_INT	owned )

Definition at line 56 of file symbfac.c.

{ css_t *css;
 
  mymalloc(css, 1, css_t);
  mymalloc(css->xnzl, (neqs+1), PORD_INT);
  mymalloc(css->xnzlsub, neqs, PORD_INT);
  if (owned)
   { mymalloc(css->nzlsub, nind, PORD_INT); }
  else
   { css->nzlsub = NULL; }
  css->neqs = neqs;
  css->nind = nind;
  css->owned = owned;
 
  return(css);
}

◆ newDenseMtx()

denseMtx_t * newDenseMtx	(	workspace_t *	,
		PORD_INT	)

◆ newDomainDecomposition()

domdec_t * newDomainDecomposition	(	PORD_INT	nvtx,
		PORD_INT	nedges )

Definition at line 100 of file ddcreate.c.

{ domdec_t *dd;
 
  mymalloc(dd, 1, domdec_t);
  mymalloc(dd->vtype, nvtx, PORD_INT);
  mymalloc(dd->color, nvtx, PORD_INT);
  mymalloc(dd->map, nvtx, PORD_INT);
 
  dd->G = newGraph(nvtx, nedges);
  dd->ndom = dd->domwght = 0;
  dd->cwght[GRAY] = dd->cwght[BLACK] = dd->cwght[WHITE] = 0;
  dd->prev = dd->next = NULL;
 
  return(dd);
}

◆ newElimGraph()

gelim_t * newElimGraph	(	PORD_INT	nvtx,
		PORD_INT	nedges )

Definition at line 108 of file gelim.c.

{ gelim_t *Gelim;
 
  mymalloc(Gelim, 1, gelim_t);
  Gelim->G = newGraph(nvtx, nedges);
  Gelim->maxedges = nedges;
 
  mymalloc(Gelim->len, nvtx, PORD_INT);
  mymalloc(Gelim->elen, nvtx, PORD_INT);
  mymalloc(Gelim->parent, nvtx, PORD_INT);
  mymalloc(Gelim->degree, nvtx, PORD_INT);
  mymalloc(Gelim->score, nvtx, PORD_INT);
 
  return(Gelim);
}

◆ newElimTree()

elimtree_t * newElimTree	(	PORD_INT	nvtx,
		PORD_INT	nfronts )

Definition at line 110 of file tree.c.

{ elimtree_t *T;
 
  mymalloc(T, 1, elimtree_t);
  mymalloc(T->ncolfactor, nfronts, PORD_INT);
  mymalloc(T->ncolupdate, nfronts, PORD_INT);
  mymalloc(T->parent, nfronts, PORD_INT);
  mymalloc(T->firstchild, nfronts, PORD_INT);
  mymalloc(T->silbings, nfronts, PORD_INT);
  mymalloc(T->vtx2front, nvtx, PORD_INT);
 
  T->nvtx = nvtx;
  T->nfronts = nfronts;
  T->root = -1;
 
  return(T);
}

◆ newFactorMtx()

factorMtx_t * newFactorMtx ( PORD_INT nelem )

Definition at line 447 of file symbfac.c.

{ factorMtx_t *L;
 
  mymalloc(L, 1, factorMtx_t);
  mymalloc(L->nzl, nelem, FLOAT);
 
  L->nelem = nelem;
  L->css = NULL;
  L->frontsub = NULL;
  L->perm = NULL;
 
  return(L);
}

◆ newFrontSubscripts()

frontsub_t * newFrontSubscripts ( elimtree_t * PTP )

Definition at line 265 of file symbfac.c.

{ frontsub_t *frontsub;
  PORD_INT        nfronts, nind;
 
  nfronts = PTP->nfronts;
  nind = nFactorIndices(PTP);
 
  mymalloc(frontsub, 1, frontsub_t);
  mymalloc(frontsub->xnzf, (nfronts+1), PORD_INT);
  mymalloc(frontsub->nzfsub, nind, PORD_INT);
 
  frontsub->PTP = PTP;
  frontsub->nind = nind;
 
  return(frontsub);
}

◆ newFrontSubscriptsPAR()

frontsub_t * newFrontSubscriptsPAR	(	mask_t *	,
		mapping_t *	,
		elimtree_t *	)

◆ newGbisect()

gbisect_t * newGbisect ( graph_t * G )

Definition at line 59 of file gbisect.c.

{ gbisect_t *Gbisect;
 
  mymalloc(Gbisect, 1, gbisect_t);
  mymalloc(Gbisect->color, G->nvtx, PORD_INT);
 
  Gbisect->G = G;
  Gbisect->cwght[GRAY] = 0;
  Gbisect->cwght[BLACK] = 0;
  Gbisect->cwght[WHITE] = 0;
 
  return(Gbisect);
}

◆ newGraph()

graph_t * newGraph	(	PORD_INT	nvtx,
		PORD_INT	nedges )

Definition at line 50 of file graph.c.

{ graph_t *G;
  PORD_INT     i;
 
  mymalloc(G, 1, graph_t);
  mymalloc(G->xadj, (nvtx+1), PORD_INT);
  mymalloc(G->adjncy, nedges, PORD_INT);
  mymalloc(G->vwght, nvtx, PORD_INT);
 
  G->nvtx = nvtx;
  G->nedges = nedges;
  G->type = UNWEIGHTED;
  G->totvwght = nvtx;
  for (i = 0; i < nvtx; i++)
    G->vwght[i] = 1;
 
  return(G);
}

◆ newInputMtx()

inputMtx_t * newInputMtx	(	PORD_INT	,
		PORD_INT	)

◆ newMapping()

mapping_t * newMapping	(	elimtree_t *	,
		PORD_INT	)

◆ newMask()

mask_t * newMask ( PORD_INT )

◆ newMinPriority()

minprior_t * newMinPriority	(	PORD_INT	nvtx,
		PORD_INT	nstages )

Definition at line 81 of file minpriority.c.

{ minprior_t  *minprior;
  stageinfo_t *stageinfo;
 
  mymalloc(stageinfo, nstages, stageinfo_t);
  mymalloc(minprior, 1, minprior_t);
  minprior->Gelim = NULL;
  minprior->ms = NULL;
  minprior->bucket = NULL;
  minprior->stageinfo = stageinfo;
 
  mymalloc(minprior->reachset, nvtx, PORD_INT);
  mymalloc(minprior->auxaux, nvtx, PORD_INT);
  mymalloc(minprior->auxbin, nvtx, PORD_INT);
  mymalloc(minprior->auxtmp, nvtx, PORD_INT);
 
  minprior->nreach = 0;
  minprior->flag = 1;
 
  return(minprior);
}

◆ newMultisector()

multisector_t * newMultisector ( graph_t * G )

Definition at line 68 of file multisector.c.

{ multisector_t *ms;
 
  mymalloc(ms, 1, multisector_t);
  mymalloc(ms->stage, G->nvtx, PORD_INT);
 
  ms->G = G;
  ms->nstages = 0;
  ms->nnodes = 0;
  ms->totmswght = 0;
 
  return(ms);
}

◆ newNDnode()

nestdiss_t * newNDnode	(	graph_t *	G,
		PORD_INT *	map,
		PORD_INT	nvint )

Definition at line 75 of file nestdiss.c.

{ nestdiss_t *nd;
 
  mymalloc(nd, 1, nestdiss_t);
  mymalloc(nd->intvertex, nvint, PORD_INT);
  mymalloc(nd->intcolor, nvint, PORD_INT);
 
  nd->G = G;
  nd->map = map;
  nd->depth = 0;
  nd->nvint = nvint;
  nd->cwght[GRAY] = nd->cwght[BLACK] = nd->cwght[WHITE] = 0;
  nd->parent = nd->childB = nd->childW = NULL;
 
  return(nd);
}

◆ newTopology()

topology_t * newTopology ( PORD_INT )

◆ nextPostorder()

PORD_INT nextPostorder	(	elimtree_t *	T,
		PORD_INT	J )

Definition at line 243 of file tree.c.

{ PORD_INT *parent, *firstchild, *silbings;
 
  parent = T->parent;
  firstchild = T->firstchild;
  silbings = T->silbings;
 
  if (silbings[J] != -1)
   { J = silbings[J];
     while (firstchild[J] != -1)
       J = firstchild[J];
   }
  else J = parent[J];
  return(J);
}

◆ nextPreorder()

PORD_INT nextPreorder	(	elimtree_t *	T,
		PORD_INT	J )

Definition at line 272 of file tree.c.

{ PORD_INT *parent, *firstchild, *silbings;
 
  parent = T->parent;
  firstchild = T->firstchild;
  silbings = T->silbings;
 
  if (firstchild[J] != -1)
    J = firstchild[J];
  else
   { while ((silbings[J] == -1) && (parent[J] != -1))
       J = parent[J];
     J = silbings[J];
   }
  return(J);
}

◆ nFactorEntries()

PORD_INT nFactorEntries ( elimtree_t * T )

Definition at line 892 of file tree.c.

{ PORD_INT *ncolfactor, *ncolupdate;
  PORD_INT ent, tri, rec, K;
 
  ncolfactor = T->ncolfactor;
  ncolupdate = T->ncolupdate;
 
  ent = 0;
  for (K = firstPostorder(T); K != -1; K = nextPostorder(T, K))
   { tri = ncolfactor[K];
     rec = ncolupdate[K];
     ent += (tri * (tri+1)) / 2;
     ent += (tri * rec);
   }
  return(ent);
}

◆ nFactorIndices()

PORD_INT nFactorIndices ( elimtree_t * T )

Definition at line 874 of file tree.c.

{ PORD_INT *ncolfactor, *ncolupdate;
  PORD_INT nfronts, ind, K;
 
  nfronts = T->nfronts;
  ncolfactor = T->ncolfactor;
  ncolupdate = T->ncolupdate;
 
  ind = 0;
  for (K = 0; K < nfronts; K++)
    ind += (ncolfactor[K] + ncolupdate[K]);
  return(ind);
}

◆ nFactorOps()

FLOAT nFactorOps ( elimtree_t * T )

Definition at line 913 of file tree.c.

{ PORD_INT   *ncolfactor, *ncolupdate;
  FLOAT ops, tri, rec;
  PORD_INT   K;
 
  ncolfactor = T->ncolfactor;
  ncolupdate = T->ncolupdate;
 
  ops = 0.0;
  for (K = firstPostorder(T); K != -1; K = nextPostorder(T, K))
   { tri = ncolfactor[K];
     rec = ncolupdate[K];
     ops += (tri*tri*tri) / 3.0 + (tri*tri) / 2.0 - (5*tri) / 6.0;
     ops += (tri*tri*rec) + (rec*(rec+1)*tri);
   }
  return(ops);
}

◆ nTriangularOps()

FLOAT nTriangularOps ( elimtree_t * T )

Definition at line 957 of file tree.c.

{ PORD_INT   *ncolfactor, *ncolupdate;
  FLOAT ops, tri, rec;
  PORD_INT   K;
 
  ncolfactor = T->ncolfactor;
  ncolupdate = T->ncolupdate;
 
  ops = 0.0;
  for (K = firstPostorder(T); K != -1; K = nextPostorder(T, K))
   { tri = ncolfactor[K];
     rec = ncolupdate[K];
     ops += (tri*tri) + 2.0*tri*rec;  /* forward ops */
     ops += (tri*tri) + 2.0*tri*rec;  /* backward ops */
   }
  return(ops);
}

◆ numfac()

void numfac	(	factorMtx_t *	L,
		timings_t *	cpus )

◆ numfacPAR()

void numfacPAR	(	topology_t *	,
		mask_t *	,
		mapping_t *	,
		factorMtx_t *	,
		PORD_INT	msglvl,
		timings_t *	)

◆ nWorkspace()

PORD_INT nWorkspace ( elimtree_t * T )

Definition at line 817 of file tree.c.

{ PORD_INT *ncolfactor, *ncolupdate, *firstchild, *silbings, *minWspace;
  PORD_INT nfronts, K, ncolfrontK, frontsizeK, wspace, child, nxtchild, m, s;
 
  nfronts = T->nfronts;
  ncolfactor = T->ncolfactor;
  ncolupdate = T->ncolupdate;
  firstchild = T->firstchild;
  silbings = T->silbings;
 
  /* -------------------------
     set up the working arrays
     ------------------------- */
  mymalloc(minWspace, nfronts, PORD_INT);
 
  /* -------------------------------------------
     postorder traversal of the elimination tree
     ------------------------------------------- */
  wspace = 0;
  for (K = firstPostorder(T); K != -1; K = nextPostorder(T, K))
   { ncolfrontK = ncolfactor[K] + ncolupdate[K];
     frontsizeK = (ncolfrontK * (ncolfrontK + 1)) >> 1;
 
     if ((child = firstchild[K]) == -1)
       minWspace[K] = frontsizeK;
     else
      { child = firstchild[K];
        nxtchild = silbings[child];
        m = s = minWspace[child];
        while (nxtchild != -1)
         { s = s - minWspace[child]
               + ((ncolupdate[child] * (ncolupdate[child] + 1)) >> 1)
               + minWspace[nxtchild];
           m = max(m, s);
           child = nxtchild;
           nxtchild = silbings[nxtchild];
         }
        s = s - minWspace[child]
            + ((ncolupdate[child] * (ncolupdate[child] + 1)) >> 1)
            + frontsizeK;
        minWspace[K] = max(m, s);
      }
 
     wspace = max(wspace, minWspace[K]);
   }
 
  /* ----------------------
     free memory and return
     ---------------------- */
  free(minWspace);
  return(wspace);
}

◆ orderMinPriority()

elimtree_t * orderMinPriority	(	minprior_t *	minprior,
		options_t *	options,
		timings_t *	cpus )

Definition at line 159 of file minpriority.c.

{ elimtree_t *T;
  PORD_INT        nvtx, nstages, istage, scoretype, ordtype;
 
  nvtx = minprior->Gelim->G->nvtx;
  nstages = minprior->ms->nstages;
 
  ordtype = options[OPTION_ORDTYPE];
  scoretype = options[OPTION_NODE_SELECTION2];
 
  /* ------------------------------
     check whether nstages is valid
     ------------------------------ */
  if ((nstages < 1) || (nstages > nvtx))
   { fprintf(stderr, "\nError in function orderMinPriority\n"
          "  no valid number of stages in multisector (#stages = %d)\n",
          nstages);
     quit();
   }
 
  if ((nstages < 2) && (ordtype != MINIMUM_PRIORITY))
   { fprintf(stderr, "\nError in function orderMinPriority\n"
          "  not enough stages in multisector (#stages = %d)\n", nstages);
     quit();
   }
 
  /* --------------------------------------------------------------
     first stage: eliminate all vertices in the remaining subgraphs
     -------------------------------------------------------------- */
  scoretype = options[OPTION_NODE_SELECTION1];
  eliminateStage(minprior, 0, scoretype, cpus);
 
  /* -------------------------------------------------------
     other stages: eliminate all vertices in the multisector
     ------------------------------------------------------- */
  switch(ordtype)
   { case MINIMUM_PRIORITY:
       break;
     case INCOMPLETE_ND:
       for (istage = 1; istage < nstages; istage++)
         eliminateStage(minprior, istage, scoretype, cpus);
       break;
     case MULTISECTION:
       eliminateStage(minprior, nstages-1, scoretype, cpus);
       break;
     default:
       fprintf(stderr, "\nError in function orderMinPriority\n"
            "  unrecognized ordering type %d\n", ordtype);
       quit();
   }
 
  /* -------------------------------------------
     print statistics for the elimination stages
     ------------------------------------------- */
  if ((ordtype != MINIMUM_PRIORITY) && (options[OPTION_MSGLVL] > 1))
    for (istage = 0; istage < nstages; istage++)
      printf("%4d. stage: #steps %6d, weight %6d, nzl %8d, ops %e\n", istage,
             minprior->stageinfo[istage].nstep,
             minprior->stageinfo[istage].welim,
             minprior->stageinfo[istage].nzf,
             minprior->stageinfo[istage].ops);
 
  /* -----------------------------------
     extract elimination tree and return
     ----------------------------------- */
  T = extractElimTree(minprior->Gelim);
  return(T);
}

◆ permFromElimTree()

void permFromElimTree	(	elimtree_t *	T,
		PORD_INT *	perm )

Definition at line 440 of file tree.c.

{ PORD_INT *vtx2front, *first, *link;
  PORD_INT nvtx, nfronts, K, u, count;
 
  nvtx = T->nvtx;
  nfronts = T->nfronts;
  vtx2front = T->vtx2front;
 
  /* -----------------------------------------------------------
     store the vertices/columns of a front in a bucket structure
     ----------------------------------------------------------- */
  mymalloc(first, nfronts, PORD_INT);
  mymalloc(link, nvtx, PORD_INT);
 
  for (K = 0; K < nfronts; K++)
    first[K] = -1;
  for (u = nvtx-1; u >= 0; u--)
   { K = vtx2front[u];
     link[u] = first[K];
     first[K] = u;
   }
 
  /* -----------------------------------------------------
     postorder traversal of the elimination tree to obtain
     the permutation vectors perm, invp
     ----------------------------------------------------- */
  count = 0;
  for (K = firstPostorder(T); K != -1; K = nextPostorder(T, K))
    for (u = first[K]; u != -1; u = link[u])
     { perm[u] = count;
       count++;
     }
 
  /* ----------------------
     free memory and return
     ---------------------- */
  free(first); free(link);
}

◆ permuteElimTree()

elimtree_t * permuteElimTree	(	elimtree_t *	T,
		PORD_INT *	perm )

Definition at line 483 of file tree.c.

{ elimtree_t *PTP;
  PORD_INT        nvtx, nfronts, J, u;
 
  nvtx = T->nvtx;
  nfronts = T->nfronts;
 
  /* --------------------------------------------------------------
     allocate space for the new elimtree object and copy front data
     the permuted tree has the same number of fronts/tree structure
     -------------------------------------------------------------- */
  PTP = newElimTree(nvtx, nfronts);
  PTP->root = T->root;
  for (J = 0; J < nfronts; J++)
   { PTP->ncolfactor[J] = T->ncolfactor[J];
     PTP->ncolupdate[J] = T->ncolupdate[J];
     PTP->parent[J] = T->parent[J];
     PTP->firstchild[J] = T->firstchild[J];
     PTP->silbings[J] = T->silbings[J];
   }
 
  /* ---------------------------------------------------------------------
     set up the new vtx2front vector; the trees only differ in this vector
     --------------------------------------------------------------------- */
  for (u = 0; u < nvtx; u++)
    PTP->vtx2front[perm[u]] = T->vtx2front[u];
 
  return(PTP);
}

◆ permuteInputMtx()

inputMtx_t * permuteInputMtx	(	inputMtx_t *	,
		PORD_INT *	)

◆ printDenseMtx()

void printDenseMtx ( denseMtx_t * )

◆ printDomainDecomposition()

void printDomainDecomposition ( domdec_t * dd )

Definition at line 133 of file ddcreate.c.

{ graph_t *G;
  PORD_INT   count, u, v, i, istart, istop;
 
  G = dd->G;
  printf("\n#nodes %d (#domains %d, weight %d), #edges %d, totvwght %d\n",
         G->nvtx, dd->ndom, dd->domwght, G->nedges >> 1, G->totvwght);
  printf("partition weights: S %d, B %d, W %d\n", dd->cwght[GRAY],
         dd->cwght[BLACK], dd->cwght[WHITE]);
  for (u = 0; u < G->nvtx; u++)
   { count = 0;
     printf("--- adjacency list of node %d (vtype %d, color %d, map %d\n",
            u, dd->vtype[u], dd->color[u], dd->map[u]);
     istart = G->xadj[u];
     istop = G->xadj[u+1];
     for (i = istart; i < istop; i++)
      { v = G->adjncy[i];
        printf("%5d (vtype %2d, color %2d)", v, dd->vtype[v], dd->color[v]);
        if ((++count % 3) == 0)
          printf("\n");
      }
     if ((count % 3) != 0)
       printf("\n");
   }
}

◆ printElimGraph()

void printElimGraph ( gelim_t * Gelim )

Definition at line 143 of file gelim.c.

{ graph_t *G;
  PORD_INT     count, u, v, i, istart;
 
  G = Gelim->G;
  for (u = 0; u < G->nvtx; u++)
   { istart = G->xadj[u];
 
     /* ---------------------------------------------------------------
        case 1: u is a principal variable
          => vwght[u]: weight of all mapped indistinguishable variables
          => degree[u]: approximate degree
        ---------------------------------------------------------------- */
     if ((Gelim->score[u] == -1) || (Gelim->score[u] >= 0))
      { printf("--- adjacency list of variable %d (weight %d, degree %d, "
               "score %d):\n", u, G->vwght[u], Gelim->degree[u],
               Gelim->score[u]);
        printf("elements:\n");
        count = 0;
        for (i = istart; i < istart + Gelim->elen[u]; i++)
         { printf("%5d", G->adjncy[i]);
           if ((++count % 16) == 0)
             printf("\n");
         }
        if ((count % 16) != 0)
          printf("\n");
        printf("variables:\n");
        count = 0;
        for (i = istart + Gelim->elen[u]; i < istart + Gelim->len[u]; i++)
         { printf("%5d", G->adjncy[i]);
           if ((++count % 16) == 0)
             printf("\n");
         }
        if ((count % 16) != 0)
          printf("\n");
      }
 
     /* ---------------------------------------------------------------
        case 2: u is nonprincipal/removed by mass elimination
        ---------------------------------------------------------------- */
     else if (Gelim->score[u] == -2)
       printf("--- variable %d is nonprincipal/removed by mass elim. "
              "(parent %d)\n", u, Gelim->parent[u]);
 
     /* -----------------------------------------------
        case 3: u is an element:
          => degree[u]: weight of boundary
        ----------------------------------------------- */
     else if (Gelim->score[u] == -3)
      { printf("--- boundary of element %d (degree %d, score %d):"
               "\n", u, Gelim->degree[u], Gelim->score[u]);
        count = 0;
        for (i = istart; i < istart + Gelim->len[u]; i++)
         { v = G->adjncy[i];
           if (G->vwght[v] > 0)
            { printf("%5d", G->adjncy[i]);
              if ((++count % 16) == 0)
                printf("\n");
            }
         }
        if ((count % 16) != 0)
          printf("\n");
      }
 
     /* --------------------------------
        case 4: u is an absorbed element
        -------------------------------- */
     else if (Gelim->score[u] == -4)
       printf("--- element %d has been absorbed (parent %d)\n", u,
              Gelim->parent[u]);
 
     /* ----------------------------------------
        none of the above cases is true => error
        ---------------------------------------- */
     else
      { fprintf(stderr, "\nError in function printElimGraph\n"
             "  node %d has invalid score %d\n", u, Gelim->score[u]);
        quit();
      }
   }
}

◆ printElimTree()

void printElimTree ( elimtree_t * T )

Definition at line 147 of file tree.c.

{ PORD_INT *ncolfactor, *ncolupdate, *parent, *firstchild, *silbings, *vtx2front;
  PORD_INT *first, *link, nvtx, nfronts, root, J, K, u, count, child;
 
  nvtx = T->nvtx;
  nfronts = T->nfronts;
  root = T->root;
  ncolfactor = T->ncolfactor;
  ncolupdate = T->ncolupdate;
  parent = T->parent;
  firstchild = T->firstchild;
  silbings = T->silbings;
  vtx2front = T->vtx2front;
 
  printf("#fronts %d, root %d\n", nfronts, root);
 
  /* -----------------------------------------------------------
     store the vertices/columns of a front in a bucket structure
     ----------------------------------------------------------- */
  mymalloc(first, nfronts, PORD_INT);
  mymalloc(link, nvtx, PORD_INT);
 
  for (J = 0; J < nfronts; J++)
    first[J] = -1;
  for (u = nvtx-1; u >= 0; u--)
   { J = vtx2front[u];
     link[u] = first[J];
     first[J] = u;
   }
 
  /* -----------------------------------------------------------
     print fronts according to a postorder traversal of the tree
     ----------------------------------------------------------- */
  for (K = firstPostorder(T); K != -1; K = nextPostorder(T, K))
   { printf("--- front %d, ncolfactor %d, ncolupdate %d, parent %d\n",
            K, ncolfactor[K], ncolupdate[K], parent[K]);
     count = 0;
     printf("children:\n");
     for (child = firstchild[K]; child != -1; child = silbings[child])
      { printf("%5d", child);
        if ((++count % 16) == 0)
          printf("\n");
      }
     if ((count % 16) != 0)
       printf("\n");
     count = 0;
     printf("vertices mapped to front:\n");
     for (u = first[K]; u != -1; u = link[u])
      { printf("%5d", u);
        if ((++count % 16) == 0)
          printf("\n");
      }
     if ((count % 16) != 0)
       printf("\n");
   }
 
  /* ----------------------
     free memory and return
     ---------------------- */
  free(first); free(link);
}

◆ printFactorMtx()

void printFactorMtx ( factorMtx_t * L )

Definition at line 478 of file symbfac.c.

{ css_t *css;
  FLOAT *nzl;
  PORD_INT   *xnzl, *nzlsub, *xnzlsub;
  PORD_INT   neqs, nelem, nind, k, ksub, i, istart, istop;
 
  nelem = L->nelem;
  nzl = L->nzl;
  css = L->css;
 
  neqs = css->neqs;
  nind = css->nind;
  xnzl = css->xnzl;
  nzlsub = css->nzlsub;
  xnzlsub = css->xnzlsub;
 
  printf("#equations %d, #elements (+diag.) %d, #indices (+diag.) %d\n",
         neqs, nelem, nind);
  for (k = 0; k < neqs; k++)
   { printf("--- column %d\n", k);
     ksub = xnzlsub[k];
     istart = xnzl[k];
     istop = xnzl[k+1];
     for (i = istart; i < istop; i++)
       printf("  row %5d, entry %e\n", nzlsub[ksub++], nzl[i]);
   }
}

◆ printFrontSubscripts()

void printFrontSubscripts ( frontsub_t * frontsub )

Definition at line 298 of file symbfac.c.

{ elimtree_t *PTP;
  PORD_INT        *xnzf, *nzfsub, *ncolfactor, *ncolupdate, *parent;
  PORD_INT        nfronts, root, K, count, i, istart, istop;
 
  PTP = frontsub->PTP;
  xnzf = frontsub->xnzf;
  nzfsub = frontsub->nzfsub;
 
  nfronts = PTP->nfronts;
  root = PTP->root;
  ncolfactor = PTP->ncolfactor;
  ncolupdate = PTP->ncolupdate;
  parent = PTP->parent;
 
  printf("#fronts %d, root %d\n", nfronts, root);
  for (K = firstPostorder(PTP); K != -1; K = nextPostorder(PTP, K))
   { printf("--- front %d, ncolfactor %d, ncolupdate %d, parent %d\n",
            K, ncolfactor[K], ncolupdate[K], parent[K]);
     count = 0;
     istart = xnzf[K];
     istop = xnzf[K+1];
     for (i = istart; i < istop; i++)
      { printf("%5d", nzfsub[i]);
        if ((++count % 16) == 0)
          printf("\n");
      }
     if ((count % 16) != 0)
       printf("\n");
   }
}

◆ printGbipart()

void printGbipart ( gbipart_t * Gbipart )

Definition at line 91 of file gbipart.c.

{ graph_t *G;
  PORD_INT     count, u, i, istart, istop;
 
  G = Gbipart->G;
  printf("\n#vertices %d (nX %d, nY %d), #edges %d, type %d, totvwght %d\n",
         G->nvtx, Gbipart->nX, Gbipart->nY, G->nedges >> 1, G->type,
         G->totvwght);
  for (u = 0; u < G->nvtx; u++)
   { count = 0;
     printf("--- adjacency list of vertex %d (weight %d):\n", u, G->vwght[u]);
     istart = G->xadj[u];
     istop = G->xadj[u+1];
     for (i = istart; i < istop; i++)
      { printf("%5d", G->adjncy[i]);
        if ((++count % 16) == 0)
          printf("\n");
      }
     if ((count % 16) != 0)
       printf("\n");
   }
}

◆ printGbisect()

void printGbisect ( gbisect_t * Gbisect )

Definition at line 87 of file gbisect.c.

{ graph_t *G;
  PORD_INT     count, u, v, i, istart, istop;
 
  G = Gbisect->G;
  printf("\n#nodes %d, #edges %d, totvwght %d\n", G->nvtx, G->nedges >> 1,
         G->totvwght);
  printf("partition weights: S %d, B %d, W %d\n", Gbisect->cwght[GRAY],
         Gbisect->cwght[BLACK], Gbisect->cwght[WHITE]);
  for (u = 0; u < G->nvtx; u++)
   { count = 0;
     printf("--- adjacency list of node %d (weight %d, color %d)\n", u,
            G->vwght[u], Gbisect->color[u]);
     istart = G->xadj[u];
     istop = G->xadj[u+1];
     for (i = istart; i < istop; i++)
      { v = G->adjncy[i];
        printf("%5d (color %2d)", v, Gbisect->color[v]);
        if ((++count % 4) == 0)
          printf("\n");
      }
     if ((count % 4) != 0)
       printf("\n");
   }
}

◆ printGraph()

void printGraph ( graph_t * G )

Definition at line 85 of file graph.c.

{ PORD_INT count, u, i, istart, istop;
 
  printf("\n#vertices %d, #edges %d, type %d, totvwght %d\n", G->nvtx,
         G->nedges >> 1, G->type, G->totvwght);
  for (u = 0; u < G->nvtx; u++)
   { count = 0;
     printf("--- adjacency list of vertex %d (weight %d):\n", u, G->vwght[u]);
     istart = G->xadj[u];
     istop = G->xadj[u+1];
     for (i = istart; i < istop; i++)
      { printf("%5d", G->adjncy[i]);
        if ((++count % 16) == 0)
          printf("\n");
      }
     if ((count % 16) != 0)
       printf("\n");
   }
}

◆ printInputMtx()

void printInputMtx ( inputMtx_t * )

◆ printMapping()

void printMapping ( mapping_t * )

◆ printTopology()

void printTopology ( topology_t * )

◆ qsortUpFloatsWithIntKeys()

void qsortUpFloatsWithIntKeys	(	PORD_INT	n,
		FLOAT *	array,
		PORD_INT *	key,
		PORD_INT *	stack )

Definition at line 140 of file sort.c.

{ register PORD_INT i, j;
  PORD_INT   t, l, m, r, p;
  FLOAT e;
 
  l = 0; r = n-1; p = 2;
  while (p > 0)
   if ((r-l) > THRES)
    { m = l + ((r-l) >> 1);
      if (key[l] > key[r])
       { swap(array[l], array[r], e); swap(key[l], key[r], t); }
      if (key[l] > key[m])
       { swap(array[l], array[m], e); swap(key[l], key[m], t); }
      if (key[r] > key[m])
       { swap(array[m], array[r], e); swap(key[m], key[r], t); }
      m = key[r]; i = l-1; j = r;
      for (;;)
       { while (key[++i] < m);
         while (key[--j] > m);
         if (i >= j) break;
         swap(array[i], array[j], e); swap(key[i], key[j], t);
       }
      swap(array[i], array[r], e); swap(key[i], key[r], t);
      if ((i-l) > (r-i))
       { stack[p++] = l;
         stack[p++] = i-1;
         l = i+1;
       }
      else
       { stack[p++] = i+1;
         stack[p++] = r;
         r = i-1;
       }
    }
   else
    { r = stack[--p];
      l = stack[--p];
    }
  if (THRES > 0) insertUpFloatsWithIntKeys(n, array, key);
}

◆ qsortUpInts()

void qsortUpInts	(	PORD_INT	n,
		PORD_INT *	array,
		PORD_INT *	stack )

Definition at line 98 of file sort.c.

{ register PORD_INT i, j;
  PORD_INT t, l, m, r, p;
 
  l = 0; r = n-1; p = 2;
  while (p > 0)
   if ((r-l) > THRES)
    { m = l + ((r-l) >> 1);
      if (array[l] > array[r]) swap(array[l], array[r], t);
      if (array[l] > array[m]) swap(array[l], array[m], t);
      if (array[r] > array[m]) swap(array[m], array[r], t);
      m = array[r]; i = l-1; j = r;
      for (;;)
       { while (array[++i] < m);
         while (array[--j] > m);
         if (i >= j) break;
         swap(array[i], array[j], t);
       }
      swap(array[i], array[r], t);
      if ((i-l) > (r-i))
       { stack[p++] = l;
         stack[p++] = i-1;
         l = i+1;
       }
      else
       { stack[p++] = i+1;
         stack[p++] = r;
         r = i-1;
       }
    }
   else
    { r = stack[--p];
      l = stack[--p];
    }
  if (THRES > 0) insertUpInts(n, array);
}

◆ randomizeGraph()

void randomizeGraph ( graph_t * G )

Definition at line 109 of file graph.c.

{ PORD_INT *xadj, *adjncy, nvtx, u, v, len, j, i, istart, istop;
 
  nvtx = G->nvtx;
  xadj = G->xadj;
  adjncy = G->adjncy;
 
  for (u = 0; u < nvtx; u++)
   { istart = xadj[u];
     istop = xadj[u+1];
     if ((len = istop - istart) > 1)
       for (i = istart; i < istop; i++)
        { j = myrandom(len);
          swap(adjncy[i], adjncy[i+j], v);
          len--;
        }
   }
}

◆ readChacoGraph()

graph_t * readChacoGraph ( char * )

◆ readHarwellBoeingMtx()

inputMtx_t * readHarwellBoeingMtx ( char * )

◆ readoutNumFacBuffer()

void readoutNumFacBuffer	(	workspace_t *	,
		buffer_t *	,
		denseMtx_t **	)

◆ readoutSymbFacBuffer()

void readoutSymbFacBuffer	(	buffer_t *	,
		frontsub_t *	,
		PORD_INT *	)

◆ readoutTriangularBuffer()

void readoutTriangularBuffer	(	buffer_t *	,
		frontsub_t *	,
		PORD_INT *	,
		FLOAT *	)

◆ recMapCube()

void recMapCube	(	topology_t *	,
		PORD_INT	,
		PORD_INT	,
		PORD_INT	,
		PORD_INT	,
		PORD_INT	,
		PORD_INT	)

◆ recvCube()

size_t recvCube	(	topology_t *	,
		void *	,
		size_t	,
		PORD_INT	)

◆ removeBucket()

void removeBucket	(	bucket_t *	bucket,
		PORD_INT	item )

Definition at line 214 of file bucket.c.

{ PORD_INT s, nextitem, lastitem;
 
  /* ----------------------------
     check whether item in bucket
     ---------------------------- */
  if (bucket->key[item] == MAX_INT)
   { fprintf(stderr, "\nError in function removeBucket\n"
          "  item %d is not in bucket\n", item);
     quit();
   }
 
  /* -----------------------
     remove item from bucket
     ----------------------- */
  nextitem = bucket->next[item];
  lastitem = bucket->last[item];
  if (nextitem != -1)
    bucket->last[nextitem] = lastitem;
  if (lastitem != -1)
    bucket->next[lastitem] = nextitem;
  else
   { s = max(0, (bucket->key[item] + bucket->offset));
     s = min(s, bucket->maxbin);
     bucket->bin[s] = nextitem;
   }
 
  /* --------------------------------------------
     decrease nobj and mark item as being removed
     -------------------------------------------- */
  bucket->nobj--;
  bucket->key[item] = MAX_INT;
}

◆ sendCube()

void sendCube	(	topology_t *	,
		void *	,
		size_t	,
		PORD_INT	)

◆ setupBipartiteGraph()

gbipart_t * setupBipartiteGraph	(	graph_t *	G,
		PORD_INT *	bipartvertex,
		PORD_INT	nX,
		PORD_INT	nY,
		PORD_INT *	vtxmap )

Definition at line 118 of file gbipart.c.

{ gbipart_t *Gbipart;
  PORD_INT       *xadj, *adjncy, *vwght, *xadjGb, *adjncyGb, *vwghtGb;
  PORD_INT       nvtx, nedgesGb, totvwght, u, x, y, i, j, jstart, jstop, ptr;
 
  nvtx = G->nvtx;
  xadj = G->xadj;
  adjncy = G->adjncy;
  vwght = G->vwght;
 
  /* ----------------------------------------------------------------
     compute number of edges and local indices of vertices in Gbipart
     ---------------------------------------------------------------- */
  nedgesGb = 0;
  for (i = 0; i < nX+nY; i++)
   { u = bipartvertex[i];
     if ((u < 0) || (u >= nvtx))
      { fprintf(stderr, "\nError in function setupBipartiteGraph\n"
             "  node %d does not belong to graph\n", u);
        quit();
      }
     jstart = xadj[u];
     jstop = xadj[u+1];
     for (j = jstart; j < jstop; j++)
       vtxmap[adjncy[j]] = -1;
     nedgesGb += (jstop - jstart);
   }
  for (i = 0; i < nX+nY; i++)
   { u = bipartvertex[i];
     vtxmap[u] = i;
   }
 
  Gbipart = newBipartiteGraph(nX, nY, nedgesGb);
  xadjGb = Gbipart->G->xadj;
  adjncyGb = Gbipart->G->adjncy;
  vwghtGb = Gbipart->G->vwght;
 
  /* ---------------------------------
     build the induced bipartite graph
     --------------------------------- */
  totvwght = 0; ptr = 0;
  for (i = 0; i < nX; i++)
   { x = bipartvertex[i];
     xadjGb[i] = ptr;
     vwghtGb[i] = vwght[x];
     totvwght += vwght[x];
     jstart = xadj[x];
     jstop = xadj[x+1];
     for (j = jstart; j < jstop; j++)
      { y = adjncy[j];
        if (vtxmap[y] >= nX)
          adjncyGb[ptr++] = vtxmap[y];
      }
   }
  for (i = nX; i < nX+nY; i++)
   { y = bipartvertex[i];
     xadjGb[i] = ptr;
     vwghtGb[i] = vwght[y];
     totvwght += vwght[y];
     jstart = xadj[y];
     jstop = xadj[y+1];
     for (j = jstart; j < jstop; j++)
      { x = adjncy[j];
        if ((vtxmap[x] >= 0) && (vtxmap[x] < nX))
          adjncyGb[ptr++] = vtxmap[x];
      }
   }
  xadjGb[nX+nY] = ptr;
  Gbipart->G->type = G->type;
  Gbipart->G->totvwght = totvwght;
  return(Gbipart);
} 

◆ setupBucket()

bucket_t * setupBucket	(	PORD_INT	maxbin,
		PORD_INT	maxitem,
		PORD_INT	offset )

Definition at line 91 of file bucket.c.

{ bucket_t *bucket;
  PORD_INT      i, u;
 
  if (offset < 0)
   { fprintf(stderr, "\nError in function setupBucket\n"
          "  offset must be >= 0\n");
     quit();
   }
 
  bucket = newBucket(maxbin, maxitem, offset);
 
  for (i = 0; i <= maxbin; i++)
    bucket->bin[i] = -1;
  for (u = 0; u <= maxitem; u++)
    { bucket->next[u] = bucket->last[u] = -1;
      bucket->key[u] = MAX_INT;
    }
 
  return(bucket);
}

◆ setupCSSFromFrontSubscripts()

css_t * setupCSSFromFrontSubscripts ( frontsub_t * frontsub )

Definition at line 221 of file symbfac.c.

{ elimtree_t *PTP;
  css_t      *css;
  PORD_INT        *xnzf, *nzfsub, *ncolfactor, *xnzl, *xnzlsub;
  PORD_INT        nind, nvtx, K, beg, knz, firstcol, col;
 
  PTP = frontsub->PTP;
  xnzf = frontsub->xnzf;
  nzfsub = frontsub->nzfsub;
  nind = frontsub->nind;
 
  nvtx = PTP->nvtx;
  ncolfactor = PTP->ncolfactor;
 
  /* -------------------------------------------------------
     allocate storage for the compressed subscript structure
     ------------------------------------------------------- */
  css = newCSS(nvtx, nind, FALSE);
  css->nzlsub = nzfsub;
 
  xnzl = css->xnzl;
  xnzlsub = css->xnzlsub;
 
  /* ---------------------------------------
     fill the compressed subscript structure
     --------------------------------------- */ 
  xnzl[0] = 0;
  for (K = firstPostorder(PTP); K != -1; K = nextPostorder(PTP, K))
   { beg = xnzf[K];
     knz = xnzf[K+1] - beg;
     firstcol = nzfsub[beg];
     for (col = firstcol; col < firstcol + ncolfactor[K]; col++)
      { xnzlsub[col] = beg++;
        xnzl[col+1] = xnzl[col] + knz--;
      }
   }
 
  return(css);
}

◆ setupCSSFromFrontSubscriptsPAR()

css_t * setupCSSFromFrontSubscriptsPAR	(	mask_t *	,
		mapping_t *	,
		frontsub_t *	)

◆ setupCSSFromGraph()

css_t * setupCSSFromGraph	(	graph_t *	G,
		PORD_INT *	perm,
		PORD_INT *	invp )

Definition at line 90 of file symbfac.c.

{ css_t *css;
  PORD_INT   *marker, *mergelink, *indices, *tmp, *xnzl, *xnzlsub, *nzlsub;
  PORD_INT   neqs, maxmem, u, v, col, mergecol, knz, mrk, beg, end;
  PORD_INT   fast, len, k, p, e, i, istart, istop;
 
  neqs = G->nvtx;
  maxmem = 2 * neqs;
 
  /* -------------------------
     set up the working arrays
     ------------------------- */
  mymalloc(marker, neqs, PORD_INT);
  mymalloc(indices, neqs, PORD_INT);
  mymalloc(mergelink, neqs, PORD_INT);
  mymalloc(tmp, neqs, PORD_INT);
  for (k = 0; k < neqs; k++)
    marker[k] = mergelink[k] = -1;
 
  /* -------------------------------------------------------
     allocate storage for the compressed subscript structure
     ------------------------------------------------------- */
  css = newCSS(neqs, maxmem, TRUE);
 
  xnzl = css->xnzl;
  nzlsub = css->nzlsub;
  xnzlsub = css->xnzlsub;
 
  /* ------------------------------------------------------------
     main loop: determine the subdiag. row indices of each column
     ------------------------------------------------------------ */
  xnzl[0] = 0;
  beg = end = 0;
  for (k = 0; k < neqs; k++)
   { indices[0] = k;
     knz = 1;
 
     if ((mergecol = mergelink[k]) != -1)        /* is k a leaf ??? */
      { mrk = marker[mergecol];
        fast = TRUE;
      }
     else
      { mrk = k;
        fast = FALSE;
      }
 
     /* --------------------------
        original columns (indices)
        -------------------------- */
     u = invp[k];
     istart = G->xadj[u];
     istop = G->xadj[u+1];
     for (i = istart; i < istop; i++)
      { v = G->adjncy[i];
        if ((col = perm[v]) > k)
         { indices[knz++] = col;
           if (marker[col] != mrk) fast = FALSE;
         }
      }
 
     /* --------------------------
        external columns (indices)
        -------------------------- */
     if ((fast) && (mergelink[mergecol] == -1))
      { xnzlsub[k] = xnzlsub[mergecol] + 1;
        knz = xnzl[mergecol+1] - xnzl[mergecol] - 1;
      }
     else
      { for (i = 0; i < knz; i++)
          marker[indices[i]] = k;
        while (mergecol != -1)
         { len = xnzl[mergecol+1] - xnzl[mergecol];
           istart = xnzlsub[mergecol];
           istop = istart + len;
           for (i = istart; i < istop; i++)
            { col = nzlsub[i];
              if ((col > k) && (marker[col] != k))
               { indices[knz++] = col;
                 marker[col] = k;
               }
            }
           mergecol = mergelink[mergecol];
         }
        qsortUpInts(knz, indices, tmp);
 
        /* ---------------------------------------------------
           store indices in nzlsub; resize nzlsub if too small
           --------------------------------------------------- */
        beg = end;
        xnzlsub[k] = beg;
        end = beg + knz;
        if (end > maxmem)
         { maxmem += neqs;
           myrealloc(nzlsub, maxmem, PORD_INT);
         }
        len = 0;
        for (i = beg; i < end; i++)
          nzlsub[i] = indices[len++];
      }
 
     /* ----------------------------
        append column k to mergelink
        ---------------------------- */
     if (knz > 1)
      { p = xnzlsub[k]+1;
        e = nzlsub[p];
        mergelink[k] = mergelink[e];
        mergelink[e] = k;
      }
     xnzl[k+1] = xnzl[k] + knz;
   }
 
  /* -----------------------------
     end of main loop: free memory
     ----------------------------- */
  free(marker); free(indices);
  free(tmp); free(mergelink);
 
  /* ------------------------------------------------------
     finalize the compressed subscript structure and return
     ------------------------------------------------------ */
  css->nind = xnzlsub[neqs-1] + 1;
  myrealloc(nzlsub, css->nind, PORD_INT);
  css->nzlsub = nzlsub;
  return(css);
}

◆ setupElimGraph()

gelim_t * setupElimGraph ( graph_t * G )

Definition at line 229 of file gelim.c.

{ gelim_t *Gelim;
  PORD_INT     *xadj, *adjncy, *vwght, *xadjGelim, *adjncyGelim, *vwghtGelim;
  PORD_INT     *len, *elen, *parent, *degree, *score;
  PORD_INT     nvtx, nedges, deg, u, i, istart, istop;
 
  nvtx = G->nvtx;
  nedges = G->nedges;
  xadj = G->xadj;
  adjncy = G->adjncy;
  vwght = G->vwght;
 
  Gelim = newElimGraph(nvtx, nedges+nvtx);
  xadjGelim = Gelim->G->xadj;
  adjncyGelim = Gelim->G->adjncy;
  vwghtGelim = Gelim->G->vwght;
  len = Gelim->len;
  elen = Gelim->elen;
  parent = Gelim->parent;
  degree = Gelim->degree;
  score = Gelim->score;
 
  /* --------------
     copy the graph
     -------------- */
  Gelim->G->type = G->type;
  Gelim->G->totvwght = G->totvwght;
  for (u = 0; u < nvtx; u++)
   { xadjGelim[u] = xadj[u];
     vwghtGelim[u] = vwght[u];
   }
  xadjGelim[nvtx] = xadj[nvtx];
  for (i = 0; i < nedges; i++)
    adjncyGelim[i] = adjncy[i];
  Gelim->G->nedges = nedges;
 
  /* ----------------------
     initialize all vectors
     ---------------------- */
  for (u = 0; u < nvtx; u++)
   { istart = xadj[u];
     istop = xadj[u+1];
     len[u] = istop - istart;
     elen[u] = 0;
     parent[u] = -1;
     deg = 0;
 
     switch(Gelim->G->type)  /* compute the external degree of u */
      { case UNWEIGHTED:
          deg = len[u];
          break;
        case WEIGHTED:
          for (i = istart; i < istop; i++)
            deg += vwght[adjncy[i]];
          break;
        default:
          fprintf(stderr, "\nError in function setupElimGraph\n"
               "  unrecognized graph type %d\n", Gelim->G->type);
      }
     degree[u] = deg;
 
     if (len[u] == 0)        /* len(u) = 0 => adjncy list of u not in use */
       xadjGelim[u] = -1;    /* mark with -1, otherwise crunchElimGraph fails */
     score[u] = -1;
   }
 
  return(Gelim);
}

◆ setupElimTree()

elimtree_t * setupElimTree	(	graph_t *	G,
		PORD_INT *	perm,
		PORD_INT *	invp )

Definition at line 293 of file tree.c.

{ elimtree_t *T;
  css_t      *css;
  PORD_INT        *xadj, *adjncy, *vwght, *ncolfactor, *ncolupdate, *parent;
  PORD_INT        *vtx2front, *realroot, *uf_father, *uf_size;
  PORD_INT        *xnzl, *nzlsub, *xnzlsub;
  PORD_INT        nvtx, front, front2, froot, f, r, u, v, i, istart, istop;
  PORD_INT        prevlen, len, h, hsub;
 
  nvtx = G->nvtx;
  xadj = G->xadj;
  adjncy = G->adjncy;
  vwght = G->vwght;
 
  /* --------------------------
     set up the working storage
     -------------------------- */
  mymalloc(realroot, nvtx, PORD_INT);
  mymalloc(uf_father, nvtx, PORD_INT);
  mymalloc(uf_size, nvtx, PORD_INT);
 
  /* ------------------------------------------------
     allocate storage for the elimination tree object
     ------------------------------------------------ */
  T = newElimTree(nvtx, nvtx);
  ncolfactor = T->ncolfactor;
  ncolupdate = T->ncolupdate;
  parent = T->parent;
  vtx2front = T->vtx2front;
 
  /* -----------------------------
     set up the parent vector of T
     ----------------------------- */
  for (front = 0; front < nvtx; front++)
   { parent[front] = -1;
     u = invp[front];           /* only vertex u belongs to this front */
     uf_father[front] = front;  /* front forms a set in union-find structure */
     uf_size[front] = 1;        /* the set consists of a single front */
     realroot[front] = front;
     froot = front;
 
     /* run through the adjacency list of u */
     istart = xadj[u];
     istop = xadj[u+1];
     for (i = istart; i < istop; i++)
      { v = adjncy[i];
        front2 = perm[v];
        if (front2 < front)
         { r = front2;
 
           while (uf_father[r] != r)  /* find root of front2 in union-find */
             r = uf_father[r];
           while (front2 != r)        /* path compression */
            { f = front2;
              front2 = uf_father[front2];
              uf_father[f] = r;
            }
 
           f = realroot[r];           /* merge union-find sets */
           if ((parent[f] == -1) && (f != front))
            { parent[f] = front;
              if (uf_size[froot] < uf_size[r])
               { uf_father[froot] = r;
                 uf_size[r] += uf_size[froot];
                 froot = r;
               }
              else
               { uf_father[r] = froot;
                 uf_size[froot] += uf_size[r];
               }
              realroot[froot] = front;
            }
         }
      }
   }
 
  /* ---------------------------------------------
     set the vectors T->firstchild and T->silbings
     --------------------------------------------- */
  initFchSilbRoot(T);
 
  /* ----------------------------------------------------------
     set the vectors T->vtx2front, T->ncolfactor, T->ncolupdate
     ---------------------------------------------------------- */
  css = setupCSSFromGraph(G, perm, invp);
  xnzl = css->xnzl;
  nzlsub = css->nzlsub;
  xnzlsub = css->xnzlsub;
 
  prevlen = 0;
  for (front = 0; front < nvtx; front++)
   { u = invp[front];
     ncolfactor[front] = vwght[u];
     ncolupdate[front] = 0;
     vtx2front[u] = front;
     len = xnzl[front+1] - xnzl[front];
     if (prevlen - 1 == len)
       ncolupdate[front] = ncolupdate[front-1] - vwght[u];
     else
      { h = xnzlsub[front] + 1;
        for (i = 1; i < len; i++)
         { hsub = nzlsub[h++];
           v = invp[hsub];
           ncolupdate[front] += vwght[v];
         }
      }
     prevlen = len;
   }
 
  /* ----------------------
     free memory and return
     ---------------------- */
  free(css);
  free(realroot); free(uf_father); free(uf_size);
  return(T);
}

◆ setupFrontalMtx()

denseMtx_t * setupFrontalMtx	(	workspace_t *	,
		factorMtx_t *	,
		PORD_INT	)

◆ setupFrontalMtxPAR()

denseMtx_t * setupFrontalMtxPAR	(	mask_t *	,
		PORD_INT	,
		workspace_t *	,
		factorMtx_t *	,
		PORD_INT	)

◆ setupFrontSubscripts()

frontsub_t * setupFrontSubscripts	(	elimtree_t *	PTP,
		inputMtx_t *	PAP )

Definition at line 334 of file symbfac.c.

{ frontsub_t *frontsub;
  PORD_INT        *ncolfactor, *ncolupdate, *firstchild, *silbings, *vtx2front;
  PORD_INT        *xnza, *nzasub, *xnzf, *nzfsub;
  PORD_INT        *marker, *tmp, *first, *indices;
  PORD_INT        nvtx, nfronts, col, firstcol, knz;
  PORD_INT        u, i, istart, istop, K, J;
 
  nvtx = PTP->nvtx;
  nfronts = PTP->nfronts;
  ncolfactor = PTP->ncolfactor;
  ncolupdate = PTP->ncolupdate;
  firstchild = PTP->firstchild;
  silbings = PTP->silbings;
  vtx2front = PTP->vtx2front;
 
  xnza = PAP->xnza;
  nzasub = PAP->nzasub;
 
  /* -------------------------
     set up the working arrays
     ------------------------- */
  mymalloc(marker, nvtx, PORD_INT);
  mymalloc(tmp, nvtx, PORD_INT);
  mymalloc(first, nfronts, PORD_INT);
  for (i = 0; i < nvtx; i++)
    marker[i] = -1;
 
  /* --------------------------------
     find the first column of a front
     -------------------------------- */
  for (u = nvtx-1; u >= 0; u--)
   { K = vtx2front[u];
     first[K] = u;
   }
 
  /* -----------------------------------------
     allocate storage for the front subscripts
     ----------------------------------------- */
  frontsub = newFrontSubscripts(PTP);
  xnzf = frontsub->xnzf;
  nzfsub = frontsub->nzfsub;
 
  knz = 0;
  for (K = 0; K < nfronts; K++)
   { xnzf[K] = knz;
     knz += (ncolfactor[K] + ncolupdate[K]);
   }
  xnzf[K] = knz;
 
  /* -------------------------------------------
     postorder traversal of the elimination tree
     ------------------------------------------- */
  for (K = firstPostorder(PTP); K != -1; K = nextPostorder(PTP, K))
   { knz = 0;
     indices = nzfsub + xnzf[K];
     firstcol = first[K];
 
     /* -------------------------------------
        internal columns (indices) of front K
        ------------------------------------- */
     for (col = firstcol; col < firstcol + ncolfactor[K]; col++)
      { indices[knz++] = col;
        marker[col] = K;
      }
 
     /* -------------------------------------
        external columns (indices) of front K
        ------------------------------------- */
     for (J = firstchild[K]; J != -1; J = silbings[J])
      { istart = xnzf[J];
        istop = xnzf[J+1];
        for (i = istart; i < istop; i++)
         { col = nzfsub[i];
           if ((col > firstcol) && (marker[col] != K))
            { marker[col] = K;
              indices[knz++] = col;
            }
         }
      }
 
     /* -------------------------------------
        original columns (indices) of front K
        ------------------------------------- */
     for (u = 0; u < ncolfactor[K]; u++)
      { istart = xnza[firstcol + u];
        istop = xnza[firstcol + u + 1];
        for (i = istart; i < istop; i++)
         { col = nzasub[i];
           if ((col > firstcol) && (marker[col] != K))
            { marker[col] = K;
              indices[knz++] = col;
            }
         }
      }
 
     /* ----------------
        sort the indices
        ---------------- */
     qsortUpInts(knz, indices, tmp);
   }
 
  /* ----------------------
     free memory and return
     ---------------------- */
  free(marker); free(tmp); free(first);
  return(frontsub);
}

◆ setupFrontSubscriptsPAR()

frontsub_t * setupFrontSubscriptsPAR	(	topology_t *	,
		mask_t *	,
		mapping_t *	,
		elimtree_t *	,
		inputMtx_t *	)

◆ setupGraphFromMtx()

graph_t * setupGraphFromMtx ( inputMtx_t * A )

Definition at line 195 of file graph.c.

{ graph_t *G;
  PORD_INT     *xnza, *nzasub, *xadj, *adjncy;
  PORD_INT     neqs, nelem, nvtx, k, h1, h2, j, i, istart, istop;
 
  neqs = A->neqs;
  nelem = A->nelem;
  xnza = A->xnza;
  nzasub = A->nzasub;
 
  /* ------------------------------------
     allocate memory for unweighted graph
     ------------------------------------ */
  G = newGraph(neqs, 2*nelem);
  nvtx = G->nvtx;
  xadj = G->xadj;
  adjncy = G->adjncy;
 
  /* -----------------------------------------
     determine the size of each adjacency list
     ----------------------------------------- */
  for (k = 0; k < neqs; k++)
    xadj[k] = xnza[k+1] - xnza[k];
  for (k = 0; k < nelem; k++)
    xadj[nzasub[k]]++;
 
  /* -------------------------------------------------------------
     determine for each vertex where its adjacency list will start
     ------------------------------------------------------------- */
  h1 = xadj[0];
  xadj[0] = 0;
  for (k = 1; k <= nvtx; k++)
   { h2 = xadj[k];
     xadj[k] = xadj[k-1] + h1;
     h1 = h2;
   }
 
  /* ------------------------
     fill the adjacency lists
     ------------------------ */
  for (k = 0; k < neqs; k++)
   { istart = xnza[k];
     istop = xnza[k+1];
     for (i = istart; i < istop; i++)
      { j = nzasub[i];
        adjncy[xadj[k]++] = j;  /* store {k,j} in adjacency list of k */
        adjncy[xadj[j]++] = k;  /* store {j,k} in adjacency list of j */
      }
   }
 
  /* --------------------------------------------
     restore startpoint of each vertex and return
     -------------------------------------------- */
  for (k = nvtx-1; k > 0; k--)
    xadj[k] = xadj[k-1];
  xadj[0] = 0;
  return(G);
}

◆ setupGridGraph()

graph_t * setupGridGraph	(	PORD_INT	dimX,
		PORD_INT	dimY,
		PORD_INT	type )

Definition at line 258 of file graph.c.

{ graph_t *G;
  PORD_INT     *xadj, *adjncy, nvtx, nedges, knz, k;
 
  /* ---------------
     initializations
     --------------- */
  G = NULL;
  knz = 0;
  nvtx = dimX * dimY;
 
  /* ---------------------------------
     create unweighted grid/mesh graph
     --------------------------------- */
  if ((type == GRID) || (type == MESH))
   { nedges = 8                           /* for edge vertices */
              + 6 * (dimX-2 + dimY-2)     /* for border vertices */
              + 4 * (dimX-2) * (dimY-2);  /* for interior vertices */
     if (type == MESH)
       nedges += 4 * (dimX-1) * (dimY-1); /* diagonals */
 
     G = newGraph(nvtx, nedges);
     xadj = G->xadj;
     adjncy = G->adjncy;
 
     for (k = 0; k < nvtx; k++)
      { xadj[k] = knz;
        if ((k+1) % dimX > 0)              /*   / k+1-dimX  (MESH) */
         { adjncy[knz++] = k+1;            /* k - k+1       (GRID) */
           if (type == MESH)               /*   \ k+1+dimX  (MESH) */
            { if (k+1+dimX < nvtx)
                adjncy[knz++] = k+1+dimX;
              if (k+1-dimX >= 0)
                adjncy[knz++] = k+1-dimX;
            }
         }
        if (k % dimX > 0)                  /* k-1-dimX \    (MESH) */
         { adjncy[knz++] = k-1;            /*      k-1 - k  (GRID) */
           if (type == MESH)               /* k-1+dimX /    (MESH) */
            { if (k-1+dimX < nvtx)
                adjncy[knz++] = k-1+dimX;
              if (k-1-dimX >= 0)
                adjncy[knz++] = k-1-dimX;
            }
         }
        if (k+dimX < nvtx)                 /* k-dimX  (GRID) */
          adjncy[knz++] = k+dimX;          /*    |           */
        if (k-dimX >= 0)                   /*    k           */
          adjncy[knz++] = k-dimX;          /*    |           */
      }                                    /* k+dimX  (GRID) */
     xadj[nvtx] = knz;
   }
 
  /* -----------------------------
     create unweighted torus graph
     ----------------------------- */
  if (type == TORUS)
   { nedges = 4 * dimX * dimY;
 
     G = newGraph(nvtx, nedges);
     xadj = G->xadj;
     adjncy = G->adjncy;
 
     for (k = 0; k < nvtx; k++)
      { xadj[k] = knz;
        if (((k+1) % dimX) == 0)           /* k -- k+1 */
          adjncy[knz++] = k+1-dimX;
        else
          adjncy[knz++] = k+1;
        if ((k % dimX) == 0)               /* k-1 -- k */
          adjncy[knz++] = k-1+dimX;
        else
          adjncy[knz++] = k-1;
        adjncy[knz++] = (k+dimX) % nvtx;               /* k-dimX */
        adjncy[knz++] =  (k+dimX*(dimY-1)) % nvtx;     /*    |   */
      }                                                /*    k   */
     xadj[nvtx] = knz;                                 /*    |   */
   }                                                   /* k+dimX */
 
  return(G);
}

◆ setupInputMtxFromGraph()

inputMtx_t * setupInputMtxFromGraph ( graph_t * )

◆ setupLaplaceMtx()

inputMtx_t * setupLaplaceMtx	(	PORD_INT	,
		PORD_INT	,
		PORD_INT	)

◆ setupMapping()

mapping_t * setupMapping	(	elimtree_t *	,
		PORD_INT	,
		PORD_INT	)

◆ setupMask()

mask_t * setupMask	(	PORD_INT	,
		PORD_INT	,
		PORD_INT	)

◆ setupMinPriority()

minprior_t * setupMinPriority ( multisector_t * ms )

Definition at line 123 of file minpriority.c.

{ minprior_t  *minprior;
  stageinfo_t *stageinfo;
  PORD_INT         *auxbin, *auxtmp;
  PORD_INT         nvtx, nstages, istage, u;
 
  nvtx = ms->G->nvtx;
  nstages = ms->nstages;
 
  minprior = newMinPriority(nvtx, nstages);
  minprior->ms = ms;
  minprior->Gelim = setupElimGraph(ms->G);
  minprior->bucket = setupBucket(nvtx, nvtx, 0);
 
  auxbin = minprior->auxbin;
  auxtmp = minprior->auxtmp;
  for (u = 0; u < nvtx; u++)
   { auxbin[u] = -1;
     auxtmp[u] = 0;
   }
 
  for (istage = 0; istage < nstages; istage++)
   { stageinfo = minprior->stageinfo + istage;
     stageinfo->nstep = 0;
     stageinfo->welim = 0;
     stageinfo->nzf = 0;
     stageinfo->ops = 0.0;
   }
 
  return(minprior);
}

◆ setupNDroot()

nestdiss_t * setupNDroot	(	graph_t *	G,
		PORD_INT *	map )

Definition at line 107 of file nestdiss.c.

{ nestdiss_t *ndroot;
  PORD_INT        *intvertex, nvtx, i;
 
  nvtx = G->nvtx;
  ndroot = newNDnode(G, map, nvtx);
  intvertex = ndroot->intvertex;
 
  for (i = 0; i < nvtx; i++)
    intvertex[i] = i;
 
  return(ndroot);
}

◆ setupNumFacBuffer()

buffer_t * setupNumFacBuffer	(	workspace_t *	,
		mask_t *	,
		PORD_INT	)

◆ setupSubgraph()

graph_t * setupSubgraph	(	graph_t *	G,
		PORD_INT *	intvertex,
		PORD_INT	nvint,
		PORD_INT *	vtxmap )

Definition at line 132 of file graph.c.

{ graph_t *Gsub;
  PORD_INT     *xadj, *adjncy, *vwght, *xadjGsub, *adjncyGsub, *vwghtGsub;
  PORD_INT     nvtx, nedgesGsub, totvwght, u, v, i, j, jstart, jstop, ptr;
 
  nvtx = G->nvtx;
  xadj = G->xadj;
  adjncy = G->adjncy;
  vwght = G->vwght;
 
  /* -------------------------------------------------------------
     compute number of edges and local indices of vertices in Gsub
     ------------------------------------------------------------- */
  nedgesGsub = 0;
  for (i = 0; i < nvint; i++)
   { u = intvertex[i];
     if ((u < 0) || (u >= nvtx))
      { fprintf(stderr, "\nError in function setupSubgraph\n"
             "  node %d does not belong to graph\n", u);
        quit();
      }
     jstart = xadj[u];
     jstop = xadj[u+1];
     for (j = jstart; j < jstop; j++)
       vtxmap[adjncy[j]] = -1;
     nedgesGsub += (jstop - jstart);
   }
  for (i = 0; i < nvint; i++)
   { u = intvertex[i];
     vtxmap[u] = i;
   }
 
  Gsub = newGraph(nvint, nedgesGsub);
  xadjGsub = Gsub->xadj;
  adjncyGsub = Gsub->adjncy;
  vwghtGsub = Gsub->vwght;
 
  /* --------------------------
     build the induced subgraph
     -------------------------- */
  totvwght = 0; ptr = 0;
  for (i = 0; i < nvint; i++)
   { u = intvertex[i];
     xadjGsub[i] = ptr;
     vwghtGsub[i] = vwght[u];
     totvwght += vwght[u];
     jstart = xadj[u];
     jstop = xadj[u+1];
     for (j = jstart; j < jstop; j++)
      { v = adjncy[j];
        if (vtxmap[v] >= 0)
          adjncyGsub[ptr++] = vtxmap[v];
      }
   }
  xadjGsub[nvint] = ptr;
  Gsub->type = G->type;
  Gsub->totvwght = totvwght;
  return(Gsub);
}

◆ setupSymbFacBuffer()

buffer_t * setupSymbFacBuffer	(	frontsub_t *	,
		PORD_INT *	)

◆ setupTopology()

topology_t * setupTopology ( void )

◆ setupTriangularBuffer()

buffer_t * setupTriangularBuffer	(	frontsub_t *	,
		PORD_INT *	,
		FLOAT *	)

◆ setupUpdateMtxFromBuffer()

denseMtx_t * setupUpdateMtxFromBuffer	(	workspace_t *	,
		FLOAT *	)

◆ setupUpdateMtxFromFrontalMtx()

denseMtx_t * setupUpdateMtxFromFrontalMtx	(	denseMtx_t *	,
		factorMtx_t *	)

◆ setupUpdateMtxFromFrontalMtxPAR()

denseMtx_t * setupUpdateMtxFromFrontalMtxPAR	(	denseMtx_t *	,
		factorMtx_t *	)

◆ shrinkDomainDecomposition()

void shrinkDomainDecomposition	(	domdec_t *	dd1,
		PORD_INT	scoretype )

Definition at line 898 of file ddcreate.c.

{ domdec_t *dd2;
  PORD_INT      *msvtxlist, *rep, *key;
  PORD_INT      nvtxdd1, nlist, u;
 
  nvtxdd1 = dd1->G->nvtx;
  mymalloc(msvtxlist, nvtxdd1, PORD_INT);
  mymalloc(rep, nvtxdd1, PORD_INT);
  mymalloc(key, nvtxdd1, PORD_INT);
 
  /* ---------------
     initializations
     --------------- */
  nlist = 0;
  for (u = 0; u < nvtxdd1; u++)
   { if (dd1->vtype[u] == 2)
       msvtxlist[nlist++] = u;
     rep[u] = u;
   }
 
  /* -------------------------------------
     compute priorities and sort multisecs
     ------------------------------------- */
  computePriorities(dd1, msvtxlist, key, scoretype);
  distributionCounting(nlist, msvtxlist, key);
 
  /* ----------------------------------------------------------
     eliminate multisecs and build coarser domain decomposition
     ---------------------------------------------------------- */
  eliminateMultisecs(dd1, msvtxlist, rep);
  findIndMultisecs(dd1, msvtxlist, rep);
  dd2 = coarserDomainDecomposition(dd1, rep);
 
  /* -----------------------------------
     append coarser domain decomposition
     ----------------------------------- */
  dd1->next = dd2;
  dd2->prev = dd1;
 
  free(msvtxlist);
  free(rep);
  free(key);
}

◆ smoothBy2Layers()

PORD_INT smoothBy2Layers	(	gbisect_t *	Gbisect,
		PORD_INT *	bipartvertex,
		PORD_INT *	pnX,
		PORD_INT	black,
		PORD_INT	white )

Definition at line 293 of file gbisect.c.

{ gbipart_t *Gbipart;
  PORD_INT       *xadj, *adjncy, *color, *cwght, *map;
  PORD_INT       *flow, *rc, *matching, *dmflag, dmwght[6];
  PORD_INT       nvtx, smoothed, nX, nX2, nY, x, y, u, i, j, jstart, jstop;
 
  nvtx = Gbisect->G->nvtx;
  xadj = Gbisect->G->xadj;
  adjncy = Gbisect->G->adjncy;
  color = Gbisect->color;
  cwght = Gbisect->cwght;
  nX = *pnX;
 
  /* ----------------------------------------------------
     map vector identifies vertices of Gbisect in Gbipart
     ---------------------------------------------------- */
  mymalloc(map, nvtx, PORD_INT);
 
  /* ----------------------------------
     construct set Y of bipartite graph
     ---------------------------------- */
  nY = 0;
  for (i = 0; i < nX; i++)
   { x = bipartvertex[i];
     jstart = xadj[x];
     jstop = xadj[x+1];
     for (j = jstart; j < jstop; j++)
      { y = adjncy[j];
        if (color[y] == black)
         { bipartvertex[nX+nY++] = y;
           color[y] = GRAY;
         }
      }
   }
  for (i = nX; i < nX+nY; i++)
   { y = bipartvertex[i];
     color[y] = black;
   }
 
  /* --------------------------------------------
     compute the Dulmage-Mendelsohn decomposition
     -------------------------------------------- */
  Gbipart = setupBipartiteGraph(Gbisect->G, bipartvertex, nX, nY, map);
 
  mymalloc(dmflag, (nX+nY), PORD_INT);
  switch(Gbipart->G->type)
   { case UNWEIGHTED:
       mymalloc(matching, (nX+nY), PORD_INT);
       maximumMatching(Gbipart, matching);
       DMviaMatching(Gbipart, matching, dmflag, dmwght);
       free(matching);
       break;
     case WEIGHTED:
       mymalloc(flow, Gbipart->G->nedges, PORD_INT);
       mymalloc(rc, (nX+nY), PORD_INT);
       maximumFlow(Gbipart, flow, rc);
       DMviaFlow(Gbipart, flow, rc, dmflag, dmwght);
       free(flow);
       free(rc);
       break;
     default:
       fprintf(stderr, "\nError in function smoothSeparator\n"
            "  unrecognized bipartite graph type %d\n", Gbipart->G->type);
       quit();
   }
 
#ifdef DEBUG
  printf("Dulmage-Mendelsohn decomp. computed\n"
         "SI %d, SX %d, SR %d, BI %d, BX %d, BR %d\n", dmwght[SI], dmwght[SX],
         dmwght[SR], dmwght[BI], dmwght[BX], dmwght[BR]);
#endif
 
  /* -----------------------------------------------------------------------
     1st TEST: try to exchange SI with BX, i.e. nodes in SI are moved from
     the separator into white (white grows), and nodes in BX are moved from
     black into the separator (black shrinks)
     ----------------------------------------------------------------------- */
  smoothed = FALSE;
  if (F(cwght[GRAY]-dmwght[SI]+dmwght[BX], cwght[black]-dmwght[BX],
        cwght[white]+dmwght[SI]) + EPS < F(cwght[GRAY], cwght[black],
                                           cwght[white]))
   { smoothed = TRUE;
 
#ifdef DEBUG
     printf("exchange SI with BX\n");
#endif
 
     cwght[white] += dmwght[SI]; cwght[GRAY] -= dmwght[SI];
     cwght[black] -= dmwght[BX]; cwght[GRAY] += dmwght[BX];
     for (i = 0; i < nX+nY; i++)
      { u = bipartvertex[i];
        if (dmflag[map[u]] == SI)
          color[u] = white;
        if (dmflag[map[u]] == BX)
          color[u] = GRAY;
      }
   }
 
  /* -----------------------------------------------------------------------
     2nd TEST: try to exchange SR with BR, i.e. nodes in SR are moved from
     the separator into white (white grows), and nodes in BR are moved from
     black into the separator (black shrinks)
     NOTE: SR is allowed to be exchanged with BR only if SI = BX = 0 or if
           SI has been exchanged with BX (Adj(SR) is a subset of BX u BR)
     ----------------------------------------------------------------------- */
  if ((F(cwght[GRAY]-dmwght[SR]+dmwght[BR], cwght[black]-dmwght[BR],
         cwght[white]+dmwght[SR]) + EPS < F(cwght[GRAY], cwght[black],
                                            cwght[white]))
      && ((smoothed) || (dmwght[SI] == 0)))
   { smoothed = TRUE;
 
#ifdef DEBUG
     printf("exchange SR with BR\n");
#endif
 
     cwght[white] += dmwght[SR]; cwght[GRAY] -= dmwght[SR];
     cwght[black] -= dmwght[BR]; cwght[GRAY] += dmwght[BR];
     for (i = 0; i < nX+nY; i++)
      { u = bipartvertex[i];
        if (dmflag[map[u]] == SR)
          color[u] = white;
        if (dmflag[map[u]] == BR)
          color[u] = GRAY;
      }
   }
 
  /* -----------------------------------------------------
     fill bipartvertex with the nodes of the new separator
     ----------------------------------------------------- */
  nX2 = 0;
  for (i = 0; i < nX+nY; i++)
   { u = bipartvertex[i];
     if (color[u] == GRAY)
       bipartvertex[nX2++] = u;
   }
  *pnX = nX2;
 
  /* -------------------------------
     free working storage and return
     ------------------------------- */
  free(map); free(dmflag);
  freeBipartiteGraph(Gbipart);
  return(smoothed);
}

◆ smoothSeparator()

void smoothSeparator	(	gbisect_t *	Gbisect,
		options_t *	options )

Definition at line 443 of file gbisect.c.

{ PORD_INT *xadj, *adjncy, *vwght, *color, *cwght, *bipartvertex;
  PORD_INT nvtx, nX, nX2, u, x, y, a, b, i, j, jstart, jstop;
 
  nvtx = Gbisect->G->nvtx;
  xadj = Gbisect->G->xadj;
  adjncy = Gbisect->G->adjncy;
  vwght = Gbisect->G->vwght;
  color = Gbisect->color;
  cwght = Gbisect->cwght;
 
  mymalloc(bipartvertex, nvtx, PORD_INT);
 
  /* ----------------------------------------------------------
     extract the separator (store its vertices in bipartvertex)
     ---------------------------------------------------------- */
  nX = 0;
  for (u = 0; u < nvtx; u++)
    if (color[u] == GRAY)
      bipartvertex[nX++] = u;
 
  do
   { /* ---------------------------------------------------------------
        minimize the separator (i.e. minimize set X of bipartite graph)
        --------------------------------------------------------------- */
     cwght[GRAY] = nX2 = 0;
     for (i = 0; i < nX; i++)
      { x = bipartvertex[i];
        a = b = FALSE;
        jstart = xadj[x];
        jstop = xadj[x+1];
        for (j = jstart; j < jstop; j++)
         { y = adjncy[j];
           if (color[y] == WHITE) a = TRUE;
           if (color[y] == BLACK) b = TRUE;
         }
        if ((a) && (!b))
         { color[x] = WHITE; cwght[WHITE] += vwght[x]; }
        else if ((!a) && (b))
         { color[x] = BLACK; cwght[BLACK] += vwght[x]; }
        else
         { bipartvertex[nX2++] = x; cwght[GRAY] += vwght[x]; }
      }
     nX = nX2;
 
#ifdef BE_CAUTIOUS
     checkSeparator(Gbisect);
#endif
 
     /* ------------------------------------------------------------------
        smooth the unweighted/weighted separator
        first pair it with the larger set; if unsuccessful try the smaller
        ------------------------------------------------------------------ */
     if (cwght[BLACK] >= cwght[WHITE])
      { a = smoothBy2Layers(Gbisect, bipartvertex, &nX, BLACK, WHITE);
        if (!a)
          a = smoothBy2Layers(Gbisect, bipartvertex, &nX, WHITE, BLACK);
      }
     else
      { a = smoothBy2Layers(Gbisect, bipartvertex, &nX, WHITE, BLACK);
        if (!a)
          a = smoothBy2Layers(Gbisect, bipartvertex, &nX, BLACK, WHITE);
      }
     if ((options[OPTION_MSGLVL] > 2) && (a))
       printf("\t separator smoothed: S %d, B %d, W %d [cost %7.2f]\n",
              cwght[GRAY], cwght[BLACK], cwght[WHITE],
              F(cwght[GRAY], cwght[BLACK], cwght[WHITE])); 
   } while (a);
     
  free(bipartvertex);
}

◆ SPACE_cleanup()

void SPACE_cleanup	(	topology_t *	,
		mask_t *	)

◆ SPACE_mapping()

mapping_t * SPACE_mapping	(	graph_t *	,
		PORD_INT *	,
		options_t *	,
		timings_t *	)

◆ SPACE_numFac()

void SPACE_numFac	(	factorMtx_t *	,
		timings_t *	)

◆ SPACE_numFacPAR()

void SPACE_numFacPAR	(	topology_t *	,
		mask_t *	,
		mapping_t *	,
		factorMtx_t *	,
		PORD_INT	msglvl,
		timings_t *	)

◆ SPACE_ordering()

elimtree_t * SPACE_ordering	(	graph_t *	G,
		options_t *	options,
		timings_t *	cpus )

Definition at line 47 of file interface.c.

{ graph_t       *Gc;
  multisector_t *ms;
  minprior_t    *minprior;
  elimtree_t    *T, *T2;
  timings_t    cpusOrd[ORD_TIME_SLOTS];
  options_t     default_options[] = { SPACE_ORDTYPE, SPACE_NODE_SELECTION1,
                       SPACE_NODE_SELECTION2, SPACE_NODE_SELECTION3,
                       SPACE_DOMAIN_SIZE, SPACE_MSGLVL };
  PORD_INT           *vtxmap, istage, totnstep, totnzf;
  FLOAT         totops;
 
  /* --------------------------------------------------
     set default options, if no other options specified
     -------------------------------------------------- */
  if (options == NULL)
    options = default_options;
 
  /* ----------------
     reset all timers
     ---------------- */
  pord_resettimer(cpusOrd[TIME_COMPRESS]);
  pord_resettimer(cpusOrd[TIME_MS]);
  pord_resettimer(cpusOrd[TIME_MULTILEVEL]);
  pord_resettimer(cpusOrd[TIME_INITDOMDEC]);
  pord_resettimer(cpusOrd[TIME_COARSEDOMDEC]);
  pord_resettimer(cpusOrd[TIME_INITSEP]);
  pord_resettimer(cpusOrd[TIME_REFINESEP]);
  pord_resettimer(cpusOrd[TIME_SMOOTH]);
  pord_resettimer(cpusOrd[TIME_BOTTOMUP]);
  pord_resettimer(cpusOrd[TIME_UPDADJNCY]);
  pord_resettimer(cpusOrd[TIME_FINDINODES]);
  pord_resettimer(cpusOrd[TIME_UPDSCORE]);
 
  /* ------------------
     compress the graph
     ------------------ */
  pord_starttimer(cpusOrd[TIME_COMPRESS]);
  mymalloc(vtxmap, G->nvtx, PORD_INT);
  Gc = compressGraph(G, vtxmap);
  pord_stoptimer(cpusOrd[TIME_COMPRESS]);
 
  if (Gc != NULL)
   { if (options[OPTION_MSGLVL] > 0)
       printf("compressed graph constructed (#nodes %d, #edges %d)\n",
              Gc->nvtx, Gc->nedges >> 1);
   }
  else
   { Gc = G;
     free(vtxmap);
     if (options[OPTION_MSGLVL] > 0)
       printf("no compressed graph constructed\n");
   }
 
  /* -------------------
     compute multisector
     ------------------- */
  
  
  pord_starttimer(cpusOrd[TIME_MS]);
  ms = constructMultisector(Gc, options, cpusOrd);
  pord_stoptimer(cpusOrd[TIME_MS]);
    
 
  if (options[OPTION_MSGLVL] > 0)
    printf("quality of multisector: #stages %d, #nodes %d, weight %d\n",
           ms->nstages, ms->nnodes, ms->totmswght);
 
  /* ---------------------------------
     compute minimum priority ordering
     --------------------------------- */
  pord_starttimer(cpusOrd[TIME_BOTTOMUP])
  minprior = setupMinPriority(ms);
  T = orderMinPriority(minprior, options, cpusOrd);
  pord_stoptimer(cpusOrd[TIME_BOTTOMUP]);
 
  if (options[OPTION_MSGLVL] > 0)
   { totnstep = totnzf = 0;
     totops = 0.0;
     for (istage = 0; istage < ms->nstages; istage++)
      { totnstep += minprior->stageinfo[istage].nstep;
        totnzf   += minprior->stageinfo[istage].nzf;
        totops   += minprior->stageinfo[istage].ops;
      }
     printf("quality of ordering: #steps %d, nzl %d, ops %e\n", totnstep,
            totnzf, totops);
   }
 
  /* -----------------------
     expand elimination tree
     ----------------------- */
  if (Gc != G)
   { T2 = expandElimTree(T, vtxmap, G->nvtx);
     freeElimTree(T);
     freeGraph(Gc);
     free(vtxmap);
   }
  else T2 = T;
 
  /* --------------------------------------------------
     pull back timing results, if vector cpus available
     -------------------------------------------------- */
  if (cpus != NULL)
   { cpus[0]  = cpusOrd[TIME_COMPRESS];
     cpus[1]  = cpusOrd[TIME_MS];
     cpus[2]  = cpusOrd[TIME_MULTILEVEL];
     cpus[3]  = cpusOrd[TIME_INITDOMDEC];
     cpus[4]  = cpusOrd[TIME_COARSEDOMDEC];
     cpus[5]  = cpusOrd[TIME_INITSEP];
     cpus[6]  = cpusOrd[TIME_REFINESEP];
     cpus[7]  = cpusOrd[TIME_SMOOTH];
     cpus[8]  = cpusOrd[TIME_BOTTOMUP];
     cpus[9]  = cpusOrd[TIME_UPDADJNCY];
     cpus[10] = cpusOrd[TIME_FINDINODES];
     cpus[11] = cpusOrd[TIME_UPDSCORE];
   }
 
  /* ----------------------
     free memory and return
     ---------------------- */
  freeMultisector(ms);
  freeMinPriority(minprior);
  return(T2);
}

◆ SPACE_setupMask()

mask_t * SPACE_setupMask	(	topology_t *	,
		PORD_INT	)

◆ SPACE_setupTopology()

topology_t * SPACE_setupTopology ( void )

◆ SPACE_solve()

void SPACE_solve	(	inputMtx_t *	,
		FLOAT *	,
		FLOAT *	,
		options_t *	,
		timings_t *	)

◆ SPACE_solveTriangular()

void SPACE_solveTriangular	(	factorMtx_t *	L,
		FLOAT *	rhs,
		FLOAT *	xvec )

◆ SPACE_solveTriangularPAR()

void SPACE_solveTriangularPAR	(	topology_t *	,
		mask_t *	,
		mapping_t *	,
		factorMtx_t *	,
		FLOAT *	,
		FLOAT *	)

◆ SPACE_solveWithPerm()

void SPACE_solveWithPerm	(	inputMtx_t *	,
		PORD_INT *	,
		FLOAT *	,
		FLOAT *	,
		options_t *	,
		timings_t *	)

◆ SPACE_solveWithPermPAR()

void SPACE_solveWithPermPAR	(	topology_t *	top,
		mask_t *	mask,
		inputMtx_t *	A,
		PORD_INT *	perm,
		FLOAT *	rhs,
		FLOAT *	xvec,
		options_t *	options,
		timings_t *	cpus )

◆ SPACE_symbFac()

factorMtx_t * SPACE_symbFac	(	elimtree_t *	,
		inputMtx_t *	)

◆ SPACE_symbFacPAR()

factorMtx_t * SPACE_symbFacPAR	(	topology_t *	,
		mask_t *	,
		mapping_t *	,
		elimtree_t *	,
		inputMtx_t *	)

◆ SPACE_transformElimTree()

elimtree_t * SPACE_transformElimTree	(	elimtree_t *	,
		PORD_INT	)

◆ split()

void split	(	mapping_t *	,
		PORD_INT	,
		PORD_INT	,
		PORD_INT	,
		PORD_INT *	,
		PORD_INT *	,
		FLOAT *	,
		PORD_INT	)

◆ splitDenseMtxColumnWise()

void splitDenseMtxColumnWise	(	denseMtx_t *	,
		mask_t *	,
		buffer_t *	,
		PORD_INT	)

◆ splitDenseMtxRowWise()

void splitDenseMtxRowWise	(	denseMtx_t *	,
		mask_t *	,
		buffer_t *	,
		PORD_INT	)

◆ splitNDnode()

void splitNDnode	(	nestdiss_t *	nd,
		options_t *	options,
		timings_t *	cpus )

Definition at line 125 of file nestdiss.c.

{ nestdiss_t *b_nd, *w_nd;
  graph_t    *Gsub;
  gbisect_t  *Gbisect;
  PORD_INT        *map, *intvertex, *intcolor, *b_intvertex, *w_intvertex;
  PORD_INT        nvint, b_nvint, w_nvint, u, i;
 
  map = nd->map;
  nvint = nd->nvint;
  intvertex = nd->intvertex;
  intcolor = nd->intcolor;
 
  /* -------------------------------------------------------------
     extract the subgraph for which a bisection has to be computed
     ------------------------------------------------------------- */
  if (nd->G->nvtx == nd->nvint)
   { Gsub = nd->G;                    /* a hack to save time and space */
     for (u = 0; u < nd->nvint; u++)  /* but do not forget the map vector */
       map[u] = u;
   }
  else
    Gsub = setupSubgraph(nd->G, intvertex, nvint, map);
  Gbisect = newGbisect(Gsub);
 
  /* ---------------------------------
     compute the bisection for Gbisect
     --------------------------------- */
  pord_starttimer(cpus[TIME_MULTILEVEL]);
  constructSeparator(Gbisect, options, cpus);
  pord_stoptimer(cpus[TIME_MULTILEVEL]);
 
  pord_starttimer(cpus[TIME_SMOOTH]);
  if (Gbisect->cwght[GRAY] > 0)
    smoothSeparator(Gbisect, options);
  pord_stoptimer(cpus[TIME_SMOOTH]);
 
  /* ----------------------------------------
     copy the bisection back to the nd object
     ---------------------------------------- */
  b_nvint = w_nvint = 0;
  nd->cwght[GRAY] = Gbisect->cwght[GRAY];
  nd->cwght[BLACK] = Gbisect->cwght[BLACK];
  nd->cwght[WHITE] = Gbisect->cwght[WHITE];
  for (i = 0; i < nvint; i++)
   { u = intvertex[i];
     intcolor[i] = Gbisect->color[map[u]];
     switch(intcolor[i])
      { case GRAY: break;
        case BLACK: b_nvint++; break;
        case WHITE: w_nvint++; break;
        default:
          fprintf(stderr, "\nError in function splitNDnode\n"
               "  node %d has unrecognized color %d\n", u, intcolor[i]);
          quit();
      }
   }
 
  /* ------------------------------------------------------
     and now split the nd object according to the bisection
     ------------------------------------------------------ */
  b_nd = newNDnode(nd->G, map, b_nvint);
  b_intvertex = b_nd->intvertex;
  w_nd = newNDnode(nd->G, map, w_nvint);
  w_intvertex = w_nd->intvertex;
 
  b_nvint = w_nvint = 0;
  for (i = 0; i < nvint; i++)
   { u = intvertex[i];
     if (intcolor[i] == BLACK) b_intvertex[b_nvint++] = u;
     if (intcolor[i] == WHITE) w_intvertex[w_nvint++] = u;
   }
  nd->childB = b_nd; b_nd->parent = nd;
  nd->childW = w_nd; w_nd->parent = nd;
  b_nd->depth = nd->depth + 1;
  w_nd->depth = nd->depth + 1;
 
  /* -----------------
     free the subgraph
     ----------------- */
  if (Gsub != nd->G)
    freeGraph(Gsub);
  freeGbisect(Gbisect);
}

◆ subtreeFactorOps()

void subtreeFactorOps	(	elimtree_t *	T,
		FLOAT *	ops )

Definition at line 935 of file tree.c.

{ PORD_INT   *ncolfactor, *ncolupdate;
  FLOAT tri, rec;
  PORD_INT   J, K;
 
  ncolfactor = T->ncolfactor;
  ncolupdate = T->ncolupdate;
 
  for (K = firstPostorder(T); K != -1; K = nextPostorder(T, K))
   { tri = ncolfactor[K];
     rec = ncolupdate[K];
     ops[K] = (tri*tri*tri) / 3.0 + (tri*tri) / 2.0 - (5*tri) / 6.0;
     ops[K] += (tri*tri*rec) + (rec*(rec+1)*tri);
     for (J = T->firstchild[K]; J != -1; J = T->silbings[J])
       ops[K] += ops[J];
   }
}

◆ trivialMultisector()

multisector_t * trivialMultisector ( graph_t * G )

Definition at line 96 of file multisector.c.

{ multisector_t *ms;
  PORD_INT           *stage, nvtx, u;
  
  /* ----------------------------------------------------------------- 
     allocate memory for the multisector object and init. stage vector
     ----------------------------------------------------------------- */
  nvtx = G->nvtx;
  ms = newMultisector(G);
  stage = ms->stage;
 
  for (u = 0; u < nvtx; u++)
    stage[u] = 0;                      /* no vertex belongs to a separator */
 
  /* -------------------------------
     finalize the multisector object
     ------------------------------- */
  ms->nstages = 1;
  ms->nnodes = 0;
  ms->totmswght = 0;
 
  return(ms);
}

◆ updateAdjncy()

void updateAdjncy	(	gelim_t *	Gelim,
		PORD_INT *	reachset,
		PORD_INT	nreach,
		PORD_INT *	tmp,
		PORD_INT *	pflag )

Definition at line 495 of file gelim.c.

{ PORD_INT *xadj, *adjncy, *vwght, *len, *elen, *parent, *score;
  PORD_INT u, v, e, me, i, j, jj, jdest, jfirstolde, jfirstv, jstart, jstop;
  PORD_INT covered, marku;
 
  xadj = Gelim->G->xadj;
  adjncy = Gelim->G->adjncy;
  vwght = Gelim->G->vwght;
  len = Gelim->len;
  elen = Gelim->elen;
  parent = Gelim->parent;
  score = Gelim->score;
 
  /* -----------------------------------------------------------------
     build the new element/variable list for each variable in reachset
     ----------------------------------------------------------------- */
  for (i = 0; i < nreach; i++)
   { u = reachset[i];
     vwght[u] = -vwght[u];        /* mark all variables in reachset */
     jstart = xadj[u];
     jstop = xadj[u] + len[u];
     jdest = jfirstolde = jstart;
 
#ifdef DEBUG
     printf("Updating adjacency list of node %d\n", u);
#endif
 
     /* --------------------------------------------------------
        scan the list of elements associated with variable u
        place newly formed elements at the beginning of the list
        -------------------------------------------------------- */
     for (j = jstart; j < jstart + elen[u]; j++)
      { e = adjncy[j];
 
#ifdef DEBUG
        printf("  >> element %d (score %d, parent %d)\n", e,score[e],parent[e]);
#endif
 
        if (score[e] == -4)       /* e has been absorbed in this elim. step */
         { me = parent[e];        /* me is the newly formed element */
           if (tmp[me] < *pflag)
            { adjncy[jdest++] = adjncy[jfirstolde];  /* move 1st old e to end */
              adjncy[jfirstolde++] = me;             /* append me at the beg. */
              tmp[me] = *pflag;
            }
         }
        else                      /* e has not been absorbed, i.e. it is */
          if (tmp[e] < *pflag)    /* an old element */
           { adjncy[jdest++] = e;
             tmp[e] = *pflag;
           }
      }
     jfirstv = jdest;             /* list of variables starts here */
 
     /* -------------------------------------------------------
        scan the list of variables associated with variable u
        place newly formed elements at the begining of the list
        ------------------------------------------------------- */
     for (j = jstart + elen[u]; j < jstop; j++)
      { v = adjncy[j];
 
#ifdef DEBUG
        printf("  >> variable %d (score %d)\n", v, score[v]);
#endif
 
        if (score[v] == -3)       /* v has been eliminated in this step */
         { if (tmp[v] < *pflag)   /* and, thus, forms a newly created elem. */
           { adjncy[jdest++] = adjncy[jfirstv];      /* move 1st var. to end  */
             adjncy[jfirstv++] = adjncy[jfirstolde]; /* move 1st old e to end */
             adjncy[jfirstolde++] = v;               /* append v at the beg.  */
             tmp[v] = *pflag;
           }
         }
        else
          adjncy[jdest++] = v;    /* v is still a variable */
      }
     elen[u] = jfirstv - jstart;
     len[u] = jdest - jstart;
     (*pflag)++;                  /* clear tmp for next round */
 
#ifdef DEBUG
     printf(" node %d: neighboring elements:\n", u);
     for (j = jstart; j < jstart + elen[u]; j++)
       printf("%5d", adjncy[j]);
     printf("\n node %d: neighboring variables:\n", u);
     for (j = jstart + elen[u]; j < jstart + len[u]; j++)
       printf("%5d", adjncy[j]);
     printf("\n");
#endif
   }
 
  /* ---------------------------------------------------------
     remove from each list all covered edges between variables
     --------------------------------------------------------- */
  for (i = 0; i < nreach; i++)
   { u = reachset[i];
     jstart = xadj[u];
     jstop = jstart + len[u];
     marku = FALSE;
 
     for (jdest = j = jstart + elen[u]; j < jstop; j++)
      { v = adjncy[j];
        if (vwght[v] > 0)         /* v does not belong to reachset */
          adjncy[jdest++] = v;    /* edge (u,v) not covered */
        if (vwght[v] < 0)         /* both vertices belong to reachset */
         { covered = FALSE;       /* check for a common element */
           if (!marku)
            { for (jj = jstart; jj < jstart + elen[u]; jj++)  /* mark elem. */
                tmp[adjncy[jj]] = *pflag;                     /* of u       */
              marku = TRUE;
            }
           for (jj = xadj[v]; jj < xadj[v] + elen[v]; jj++)   /* check elem. */
             if (tmp[adjncy[jj]] == *pflag)                   /* of v        */
              { covered = TRUE;
                break;
              }
           if (!covered)
             adjncy[jdest++] = v;
         }
      }
     len[u] = jdest - jstart;
     (*pflag)++;                  /* clear tmp for next round */
 
#ifdef DEBUG
     printf(" node %d: neighboring uncovered variables:\n", u);
     for (j = jstart + elen[u]; j < jstart + len[u]; j++)
       printf("%5d", adjncy[j]);
     printf("\n");
#endif
   }
 
  /* --------------------------------
     unmark all variables in reachset
     -------------------------------- */
  for (i = 0; i < nreach; i++)
   { u = reachset[i]; 
     vwght[u] = -vwght[u];
   }
}

◆ updateB2W()

void updateB2W	(	bucket_t *	w_bucket,
		bucket_t *	b_bucket,
		domdec_t *	dd,
		PORD_INT	domain,
		PORD_INT *	tmp_color,
		PORD_INT *	deltaW,
		PORD_INT *	deltaB,
		PORD_INT *	deltaS )

Definition at line 383 of file ddbisect.c.

{ PORD_INT *xadj, *adjncy, *vwght, *vtype;
  PORD_INT weight, u, v, i, istart, istop, j, jstart, jstop;
 
  xadj = dd->G->xadj;
  adjncy = dd->G->adjncy;
  vwght = dd->G->vwght;
  vtype = dd->vtype;
 
  istart = xadj[domain];
  istop = xadj[domain+1];
  for (i = istart; i < istop; i++)
   { u = adjncy[i];
     weight = vwght[u];
     jstart = xadj[u];
     jstop = xadj[u+1];
 
     /* ---------------------------------------------------------------
        subcase (1): before flipping domain to WHITE there was only one
          other WHITE domain v. update deltaB[v] and deltaS[v]
        --------------------------------------------------------------- */
     if (deltaW[u] < 0)
      { v = -(deltaW[u]+1);
        deltaW[u] = 1;
 
#ifdef DEBUG
        printf(" B2W case (1): (via multisec %d) removing domain %d from "
               "w_bucket\n", u, v);
#endif
 
        removeBucket(w_bucket, v);
        deltaB[v] -= weight; deltaS[v] += weight;
        insertBucket(w_bucket, deltaS[v], v);
      }
 
     /* ---------------------------------------------------------------
        subcase (2): all other domains are BLACK. update deltaB, deltaS
          of these BLACK domains. NOTE: subcase (3) may directly follow
        --------------------------------------------------------------- */
     if (deltaW[u] == 0)
      { tmp_color[u] = GRAY;
        for (j = jstart; j < jstop; j++)
         { v = adjncy[j];
           if (vtype[v] == 1)
            {
#ifdef DEBUG
              printf(" B2W case (2): (via multisec %d) removing domain %d from "
                     "b_bucket\n", u, v);
#endif
 
              removeBucket(b_bucket, v);
              deltaB[v] += weight; deltaS[v] -= weight;
              insertBucket(b_bucket, deltaS[v], v);
            }
         }
      }
 
     if (deltaB[u] < 0) deltaB[u] = 1;   /* the unique BLACK dom. flipped */
     deltaB[u]--; deltaW[u]++;
 
     /* -------------------------------------------------------------
        subcase (3): after flipping domain to WHITE there is only one
          remaining BLACK domain. search it and update deltaW, deltaS
          furthermore, store the remaining BLACK domain in deltaB[u]
        ------------------------------------------------------------- */
     if (deltaB[u] == 1)
      { for (j = jstart; j < jstop; j++)
         { v = adjncy[j];
           if ((tmp_color[v] == BLACK) && (vtype[v] == 1))
            {
#ifdef DEBUG
              printf(" B2W case (3): (via multisec %d) removing domain %d from "
                     "b_bucket\n", u, v);
#endif
 
              removeBucket(b_bucket, v);
              deltaW[v] += weight; deltaS[v] -= weight;
              deltaB[u] = -(v+1);
              insertBucket(b_bucket, deltaS[v], v);
            }
         }
      }
 
     /* -------------------------------------------------------------
        subcase (4): after flipping domain to WHITE there is no other
          BLACK domain. update deltaW, deltaS of the WHITE domains
        ------------------------------------------------------------- */
     if (deltaB[u] == 0)
      { tmp_color[u] = WHITE;
        for (j = jstart; j < jstop; j++)
         { v = adjncy[j];
           if (vtype[v] == 1)
            {
#ifdef DEBUG
              printf(" B2W case (4): (via multisec %d) removing domain %d from "
                     "w_bucket\n", u, v);
#endif
 
              removeBucket(w_bucket, v);
              deltaW[v] -= weight; deltaS[v] += weight;
              insertBucket(w_bucket, deltaS[v], v);
            }
         }
      }
   }
}

◆ updateDegree()

void updateDegree	(	gelim_t *	Gelim,
		PORD_INT *	reachset,
		PORD_INT	nreach,
		PORD_INT *	bin )

Definition at line 761 of file gelim.c.

{ PORD_INT *xadj, *adjncy, *vwght, *len, *elen, *degree;
  PORD_INT totvwght, deg, vwghtv, u, v, w, e, me, r, i, istart, istop;
  PORD_INT j, jstart, jstop;
 
  totvwght = Gelim->G->totvwght;
  xadj = Gelim->G->xadj;
  adjncy = Gelim->G->adjncy;
  vwght = Gelim->G->vwght;
  len = Gelim->len;
  elen = Gelim->elen;
  degree = Gelim->degree;
 
  /* -------------------------------------------------------------------
     degree update only for those vertices in reachset that are adjacent
     to an element
     ------------------------------------------------------------------- */
  for (r = 0; r < nreach; r++)
   { u = reachset[r];
     if (elen[u] > 0)
       bin[u] = 1;
   }
 
  /* -----------------------------------------
     and now do the approximate degree updates
     ----------------------------------------- */
  for (r = 0; r < nreach; r++)
   { u = reachset[r];
     if (bin[u] == 1)          /* me is the most recently formed element */
      { me = adjncy[xadj[u]];  /* in the neighborhood of u */
 
#ifdef DEBUG
        printf("Updating degree of all variables in L(%d) (initiated by %d)\n",
               me, u);
#endif
 
        /* ----------------------------------------------------------------
           compute in bin[e] the size of Le\Lme for all unabsorbed elements
           ---------------------------------------------------------------- */
        istart = xadj[me];
        istop = istart + len[me];         /* compute in bin[e] the size */
        for (i = istart; i < istop; i++)  /* of Le/Lme for all elements */
         { v = adjncy[i];                 /* e != me that are adjacent  */
           vwghtv = vwght[v];             /* to a principal var. e Lme  */
           if (vwghtv > 0)
            { jstart = xadj[v];
              jstop = jstart + elen[v];
              for (j = jstart; j < jstop; j++)
               { e = adjncy[j];
                 if (e != me)
                  { if (bin[e] > 0) bin[e] -= vwghtv;
                    else bin[e] = degree[e] - vwghtv;
                  }
               }
            }
         }
 
#ifdef DEBUG
        for (i = istart; i < istop; i++)
         { v = adjncy[i];
           if (vwght[v] > 0)
             for (j = xadj[v]; j < xadj[v] + elen[v]; j++)
              { e = adjncy[j];
                if (e != me)
                  printf("  >> element %d: degree %d, outer degree %d\n", e,
                         degree[e], bin[e]);
              }
         }
#endif
 
        /* ------------------------------------------------------
           update approx. degree for all v in Lme with bin[v] = 1
           ------------------------------------------------------ */
        for (i = istart; i < istop; i++)
         { v = adjncy[i];                  /* update the upper bound deg. */
           vwghtv = vwght[v];              /* of all principal variables  */
           deg = 0;                        /* in Lme that have not been   */
           if (bin[v] == 1)                /* updated yet                 */
            { jstart = xadj[v];
              jstop = jstart + len[v];
 
              /* scan the element list associated with principal v */
              for (j = jstart; j < jstart + elen[v]; j++)
               { e = adjncy[j];            
                 if (e != me) deg += bin[e];
               }
 
              /* scan the supervariables in the list associated with v */
              for (j = jstart + elen[v]; j < jstop; j++)
               { w = adjncy[j];
                 deg += vwght[w];
               }
 
              /* compute the external degree of v (add size of Lme) */
              deg = min(degree[v], deg);
              degree[v] = max(1, min(deg+degree[me]-vwghtv, totvwght-vwghtv));
              bin[v] = -1;
 
#ifdef DEBUG
              printf("  >> variable %d (totvwght %d, vwght %d): deg %d, "
                     "degme %d, approx degree %d\n", v, totvwght, vwghtv, deg,
                     degree[me], degree[v]);
#endif
            }
         }
 
        /* ------------------------------------
           clear bin[e] of all elements e != me
           ------------------------------------ */
        for (i = istart; i < istop; i++)
         { v = adjncy[i];
           vwghtv = vwght[v];
           if (vwghtv > 0)
            { jstart = xadj[v];
              jstop = jstart + elen[v];
              for (j = jstart; j < jstop; j++)
               { e = adjncy[j];
                 if (e != me) bin[e] = -1;
               }
            }
         }
      }
   }
}

◆ updateScore()

void updateScore	(	gelim_t *	Gelim,
		PORD_INT *	reachset,
		PORD_INT	nreach,
		PORD_INT	scoretype,
		PORD_INT *	bin )

Definition at line 890 of file gelim.c.

{ PORD_INT *xadj, *adjncy, *vwght, *len, *elen, *degree, *score;
  PORD_INT vwghtv, deg, degme, u, v, me, r, i, istart, istop;
  /* Modified by JYL, 16 march 2005.
   * scr could overflow for quasi dense rows.
   * Use a double instead for large degrees
   * aset it near to MAX_INT in case of problem.
   */
  double scr_dbl;
  PORD_INT scr;
 
  xadj = Gelim->G->xadj;
  adjncy = Gelim->G->adjncy;
  vwght = Gelim->G->vwght;
  len = Gelim->len;
  elen = Gelim->elen;
  degree = Gelim->degree;
  score = Gelim->score;
 
  /* ------------------------------------------------------------------
     score update only for those vertices in reachset that are adjacent
     to an element
     ------------------------------------------------------------------ */
  for (r = 0; r < nreach; r++)
   { u = reachset[r];
     if (elen[u] > 0)
       bin[u] = 1;
   }
 
  /* ----------------------------
     and now do the score updates
     ---------------------------- */
  scoretype = scoretype % 10;
  for (r = 0; r < nreach; r++)
   { u = reachset[r];
     if (bin[u] == 1)          /* me is the most recently formed element */
      { me = adjncy[xadj[u]];  /* in the neighborhood of u */
 
#ifdef DEBUG
        printf("Updating score of all variables in L(%d) (initiated by %d)\n",
               me, u);
#endif
 
        istart = xadj[me];
        istop = xadj[me] + len[me];
        for (i = istart; i < istop; i++)
         { v = adjncy[i];                  /* update score of all principal */
           if (bin[v] == 1)                /* variables in Lme that have not */
            { vwghtv = vwght[v];           /* been updated yet */
              deg = degree[v];
              degme = degree[me] - vwghtv;
          if (deg > 40000 || degme > 40000)
          {
              switch(scoretype)
               { case AMD:
                   scr_dbl = (double)deg;
                   break;
                 case AMF:
                   scr_dbl = (double)deg*(double)(deg-1)/2 - (double)degme*(double)(degme-1)/2;
                   break;
                 case AMMF:
                   scr_dbl = ((double)deg*(double)(deg-1)/2 - (double)degme*(double)(degme-1)/2) / (double)vwghtv;
                   break;
                 case AMIND:
                   scr_dbl = max(0, ((double)deg*(double)(deg-1)/2 - (double)degme*(double)(degme-1)/2)
                             - (double)deg*(double)vwghtv);
                   break;
                 default:
                   fprintf(stderr, "\nError in function updateScore\n"
                        "  unrecognized selection strategy %d\n", scoretype);
                   quit();
               }
              /* Some buckets have offset nvtx / 2.
           * Using MAX_INT - nvtx should then be safe */
              score[v] = (PORD_INT) (min(scr_dbl,MAX_INT-Gelim->G->nvtx));
          }
          else
          {
              switch(scoretype)
               { case AMD:
                   scr = deg;
                   break;
                 case AMF:
                   scr = deg*(deg-1)/2 - degme*(degme-1)/2;
                   break;
                 case AMMF:
                   scr = (deg*(deg-1)/2 - degme*(degme-1)/2) / vwghtv;
                   break;
                 case AMIND:
                   scr = max(0, (deg*(deg-1)/2 - degme*(degme-1)/2)
                             - deg*vwghtv);
                   break;
                 default:
                   fprintf(stderr, "\nError in function updateScore\n"
                        "  unrecognized selection strategy %d\n", scoretype);
                   quit();
                 }
               score[v] = scr;
          }
              bin[v] = -1;
 
#ifdef DEBUG
              printf("  >> variable %d (me %d): weight %d, (ext)degme %d, "
                     "degree %d, score %d\n", u, me, vwghtv, degme, degree[v],
                     score[v]);
#endif
         
              if (score[v] < 0)
               { fprintf(stderr, "\nError in function updateScore\n"
                      " score[%d] = %d is negative\n", v, score[v]);
                 quit();
               }
            }
         }
      }
   }
}

◆ updateW2B()

void updateW2B	(	bucket_t *	w_bucket,
		bucket_t *	b_bucket,
		domdec_t *	dd,
		PORD_INT	domain,
		PORD_INT *	tmp_color,
		PORD_INT *	deltaW,
		PORD_INT *	deltaB,
		PORD_INT *	deltaS )

Definition at line 495 of file ddbisect.c.

{ PORD_INT *xadj, *adjncy, *vwght, *vtype;
  PORD_INT weight, u, v, i, istart, istop, j, jstart, jstop;
 
  xadj = dd->G->xadj;
  adjncy = dd->G->adjncy;
  vwght = dd->G->vwght;
  vtype = dd->vtype;
 
  istart = xadj[domain];
  istop = xadj[domain+1];
  for (i = istart; i < istop; i++)
   { u = adjncy[i];
     weight = vwght[u];
     jstart = xadj[u];
     jstop = xadj[u+1];
 
     /* ---------------------------------------------------------------
        subcase (1): before flipping domain to BLACK there was only one
          other BLACK domain v. update deltaW[v] and deltaS[v]
        --------------------------------------------------------------- */
     if (deltaB[u] < 0)
      { v = -(deltaB[u]+1);
        deltaB[u] = 1;
 
#ifdef DEBUG
        printf(" W2B case (1): (via multisec %d) removing domain %d from "
               "b_bucket\n", u, v);
#endif
 
        removeBucket(b_bucket, v);
        deltaW[v] -= weight; deltaS[v] += weight;
        insertBucket(b_bucket, deltaS[v], v);
      }
 
     /* ---------------------------------------------------------------
        subcase (2): all other domains are WHITE. update deltaW, deltaS
          of these WHITE domains. NOTE: subcase (3) may directly follow
        --------------------------------------------------------------- */
     if (deltaB[u] == 0)
      { tmp_color[u] = GRAY;
        for (j = jstart; j < jstop; j++)
         { v = adjncy[j];
           if (vtype[v] == 1)
            {
#ifdef DEBUG
              printf(" W2B case (2): (via multisec %d) removing domain %d from "
                     "w_bucket\n", u, v);
#endif
 
              removeBucket(w_bucket, v);
              deltaW[v] += weight; deltaS[v] -= weight;
              insertBucket(w_bucket, deltaS[v], v);
            }
         }
      }
 
     if (deltaW[u] < 0) deltaW[u] = 1;   /* the unique WHITE dom. flipped */
     deltaB[u]++; deltaW[u]--;
 
     /* -------------------------------------------------------------
        subcase (3): after flipping domain to BLACK there is only one
          remaining WHITE domain. search it and update deltaB, deltaS
          furthermore, store the remaining WHITE domain in deltaW[u]
        ------------------------------------------------------------- */
     if (deltaW[u] == 1)
      { for (j = jstart; j < jstop; j++)
         { v = adjncy[j];
           if ((tmp_color[v] == WHITE) && (vtype[v] == 1))
            {
#ifdef DEBUG
              printf(" W2B case (3): (via multisec %d) removing domain %d from "
                     "w_bucket\n", u, v);
#endif
 
              removeBucket(w_bucket, v);
              deltaB[v] += weight; deltaS[v] -= weight;
              deltaW[u] = -(v+1);
              insertBucket(w_bucket, deltaS[v], v);
            }
         }
      }
 
     /* ---------------------------------------------------------------
        subcase (4): after flipping domain to BLACK there is no other
          WHITE domain. update deltaB, deltaS of the BLACK domains
        --------------------------------------------------------------- */
     if (deltaW[u] == 0)
      { tmp_color[u] = BLACK;
        for (j = jstart; j < jstop; j++)
         { v = adjncy[j];
           if (vtype[v] == 1)
            {
#ifdef DEBUG
              printf(" W2B case (4): (via multisec %d) removing domain %d from "
                     "b_bucket\n", u, v);
#endif
 
              removeBucket(b_bucket, v);
              deltaB[v] -= weight; deltaS[v] += weight;
              insertBucket(b_bucket, deltaS[v], v);
            }
         }
      }
   }
}

OpenRadioss 2025.1.11 OpenRadioss project

Functions