Update approx. minimum ordering method to push and keep structural empty diagonal elements to the bottom-right part of the matrix

This commit is contained in:
Gael Guennebaud 2015-03-20 16:33:48 +01:00
parent 8d9bfb3a7b
commit 4e2b18d909

View File

@ -137,22 +137,27 @@ void minimum_degree_ordering(SparseMatrix<Scalar,ColMajor,StorageIndex>& C, Perm
degree[i] = len[i]; // degree of node i degree[i] = len[i]; // degree of node i
} }
mark = internal::cs_wclear<StorageIndex>(0, 0, w, n); /* clear w */ mark = internal::cs_wclear<StorageIndex>(0, 0, w, n); /* clear w */
elen[n] = -2; /* n is a dead element */
Cp[n] = -1; /* n is a root of assembly tree */
w[n] = 0; /* n is a dead element */
/* --- Initialize degree lists ------------------------------------------ */ /* --- Initialize degree lists ------------------------------------------ */
for(i = 0; i < n; i++) for(i = 0; i < n; i++)
{ {
bool has_diag = false;
for(p = Cp[i]; p<Cp[i+1]; ++p)
if(Ci[p]==i)
{
has_diag = true;
break;
}
d = degree[i]; d = degree[i];
if(d == 0) /* node i is empty */ if(d == 1) /* node i is empty */
{ {
elen[i] = -2; /* element i is dead */ elen[i] = -2; /* element i is dead */
nel++; nel++;
Cp[i] = -1; /* i is a root of assembly tree */ Cp[i] = -1; /* i is a root of assembly tree */
w[i] = 0; w[i] = 0;
} }
else if(d > dense) /* node i is dense */ else if(d > dense || !has_diag) /* node i is dense or has no structural diagonal element */
{ {
nv[i] = 0; /* absorb i into element n */ nv[i] = 0; /* absorb i into element n */
elen[i] = -1; /* node i is dead */ elen[i] = -1; /* node i is dead */
@ -168,6 +173,10 @@ void minimum_degree_ordering(SparseMatrix<Scalar,ColMajor,StorageIndex>& C, Perm
} }
} }
elen[n] = -2; /* n is a dead element */
Cp[n] = -1; /* n is a root of assembly tree */
w[n] = 0; /* n is a dead element */
while (nel < n) /* while (selecting pivots) do */ while (nel < n) /* while (selecting pivots) do */
{ {
/* --- Select node of minimum approximate degree -------------------- */ /* --- Select node of minimum approximate degree -------------------- */