From 0a499aa2b00a04a28576e4b87947eea504525e88 Mon Sep 17 00:00:00 2001
From: Mathieu Faverge <mathieu.faverge@inria.fr>
Date: Thu, 31 Aug 2017 15:17:13 +0200
Subject: [PATCH] Add subroutines to introduce diagonal elements in graph
structure for fake factorizations
---
spm_gen_fake_values.c | 138 ++++++++++++++++++++++++++++++++++++++++--
1 file changed, 134 insertions(+), 4 deletions(-)
diff --git a/spm_gen_fake_values.c b/spm_gen_fake_values.c
index 7869c741..3c641618 100644
--- a/spm_gen_fake_values.c
+++ b/spm_gen_fake_values.c
@@ -34,8 +34,12 @@
* Array of size spm->n allocated on entry. On exit, contains the
* degree of each vertex in the spm matrix.
*
+ *******************************************************************************
+ *
+ * @return the number of diagonal elements found during the computation.
+ *
*******************************************************************************/
-static inline void
+static inline pastix_int_t
spm_compute_degrees( const pastix_spm_t *spm,
pastix_int_t *degrees )
{
@@ -43,6 +47,7 @@ spm_compute_degrees( const pastix_spm_t *spm,
pastix_int_t *colptr = spm->colptr;
pastix_int_t *rowptr = spm->rowptr;
pastix_int_t baseval;
+ pastix_int_t diagval = 0;
baseval = spmFindBase( spm );
memset( degrees, 0, spm->n * sizeof(pastix_int_t) );
@@ -65,6 +70,9 @@ spm_compute_degrees( const pastix_spm_t *spm,
degrees[i] += 1;
}
}
+ else {
+ diagval++;
+ }
}
}
break;
@@ -81,8 +89,125 @@ spm_compute_degrees( const pastix_spm_t *spm,
degrees[i] += 1;
}
}
+ else {
+ diagval++;
+ }
}
}
+
+ return diagval;
+}
+
+/**
+ *******************************************************************************
+ *
+ * @ingroup spm_dev_driver
+ *
+ * @brief Insert diagonal elements to the graph to have a full Laplacian
+ * generated
+ *
+ *******************************************************************************
+ *
+ * @param[in] spm
+ * At start, the initial spm structure with missing diagonal elements.
+ * At exit, contains the same sparse matrix with diagonal elements added.
+ *
+ * @param[in] diagval
+ * The number of diagonal elements already present in the matrix.
+ *
+ *******************************************************************************/
+static inline void
+spm_add_diag( pastix_spm_t *spm,
+ pastix_int_t diagval )
+{
+ pastix_spm_t oldspm;
+ pastix_int_t i, j, k;
+ pastix_int_t *oldcol = spm->colptr;
+ pastix_int_t *oldrow = spm->rowptr;
+ pastix_int_t *newrow, *newcol;
+ pastix_int_t baseval;
+
+ baseval = spmFindBase( spm );
+
+ memcpy( &oldspm, spm, sizeof(pastix_spm_t));
+
+ spm->nnz = oldspm.nnz + (spm->n - diagval);
+ newrow = malloc( spm->nnz * sizeof(pastix_int_t) );
+
+ switch(spm->fmttype)
+ {
+ case PastixCSR:
+ /* Swap pointers to call CSC */
+ oldcol = spm->rowptr;
+ oldrow = spm->colptr;
+ spm->colptr = newrow;
+
+ case PastixCSC:
+ newcol = oldcol;
+ if ( spm->fmttype == PastixCSC ) {
+ spm->rowptr = newrow;
+ }
+ diagval = 0;
+ for(j=0; j<spm->n; j++, newcol++) {
+ pastix_int_t nbelt = newcol[1] - newcol[0];
+ int diag = 0;
+
+ memcpy( newrow, oldrow, nbelt * sizeof(pastix_int_t) );
+ newrow += nbelt;
+
+ for(k=0; k<nbelt; k++, oldrow++) {
+ i = *oldrow - baseval;
+
+ if ( i == j ) {
+ diag = 1;
+ }
+ }
+ newcol[0] += diagval;
+ if ( !diag ) {
+ *newrow = j + baseval;
+ newrow++;
+ diagval++;
+ }
+ }
+ newcol[0] += diagval;
+
+ if ( spm->fmttype == PastixCSC ) {
+ free( oldspm.rowptr );
+ }
+ else {
+ free( oldspm.colptr );
+ }
+ assert( diagval == spm->n );
+ break;
+
+ case PastixIJV:
+ newcol = malloc( spm->nnz * sizeof(pastix_int_t) );
+ spm->colptr = newcol;
+ spm->rowptr = newrow;
+
+ for(k=0; k<spm->n; k++, newrow++, newcol++)
+ {
+ *newrow = k + baseval;
+ *newcol = k + baseval;
+ }
+
+ for(k=0; k<spm->nnz; k++, oldrow++, oldcol++)
+ {
+ i = *oldrow - baseval;
+ j = *oldcol - baseval;
+
+ if ( i == j ) {
+ continue;
+ }
+
+ *newrow = i + baseval;
+ *newcol = j + baseval;
+ newrow++;
+ newcol++;
+ }
+ free( oldspm.colptr );
+ free( oldspm.rowptr );
+ }
}
/**
@@ -178,9 +303,9 @@ spm_generate_fake_values( pastix_spm_t *spm,
void
spmGenFakeValues( pastix_spm_t *spm )
{
- pastix_int_t *degrees;
+ pastix_int_t *degrees, diagval;
double alpha = 10.;
- double beta = 1.;
+ double beta = 1.;
assert( spm->flttype == PastixPattern );
assert( spm->values == NULL );
@@ -216,7 +341,12 @@ spmGenFakeValues( pastix_spm_t *spm )
}
degrees = malloc( spm->n * sizeof(pastix_int_t));
- spm_compute_degrees( spm, degrees );
+ diagval = spm_compute_degrees( spm, degrees );
+ if ( diagval != spm->n ) {
+ /* Diagonal elements must be added to the sparse matrix */
+ spm_add_diag( spm, diagval );
+ spmUpdateComputedFields( spm );
+ }
spm_generate_fake_values( spm, degrees, alpha, beta );
free( degrees );
--
GitLab