|
kusano |
7d535a |
#include <stdio.h></stdio.h>
|
|
kusano |
7d535a |
#include "slu_ddefs.h"
|
|
kusano |
7d535a |
#include "slu_util.h"
|
|
kusano |
7d535a |
|
|
kusano |
7d535a |
|
|
kusano |
7d535a |
void
|
|
kusano |
7d535a |
dreadtriple(int *m, int *n, int *nonz,
|
|
kusano |
7d535a |
double **nzval, int **rowind, int **colptr)
|
|
kusano |
7d535a |
{
|
|
kusano |
7d535a |
/*
|
|
kusano |
7d535a |
* Output parameters
|
|
kusano |
7d535a |
* =================
|
|
kusano |
7d535a |
* (a,asub,xa): asub[*] contains the row subscripts of nonzeros
|
|
kusano |
7d535a |
* in columns of matrix A; a[*] the numerical values;
|
|
kusano |
7d535a |
* row i of A is given by a[k],k=xa[i],...,xa[i+1]-1.
|
|
kusano |
7d535a |
*
|
|
kusano |
7d535a |
*/
|
|
kusano |
7d535a |
int i, j, k, jsize, lasta, nnz, nz;
|
|
kusano |
7d535a |
double *a, *val;
|
|
kusano |
7d535a |
int *asub, *xa, *row, *col;
|
|
kusano |
7d535a |
|
|
kusano |
7d535a |
/* Matrix format:
|
|
kusano |
7d535a |
* First line: #rows, #cols, #non-zero
|
|
kusano |
7d535a |
* Triplet in the rest of lines:
|
|
kusano |
7d535a |
* row, col, value
|
|
kusano |
7d535a |
*/
|
|
kusano |
7d535a |
|
|
kusano |
7d535a |
scanf("%d%d", n, nonz);
|
|
kusano |
7d535a |
*m = *n;
|
|
kusano |
7d535a |
printf("m %d, n %d, nonz %d\n", *m, *n, *nonz);
|
|
kusano |
7d535a |
dallocateA(*n, *nonz, nzval, rowind, colptr); /* Allocate storage */
|
|
kusano |
7d535a |
a = *nzval;
|
|
kusano |
7d535a |
asub = *rowind;
|
|
kusano |
7d535a |
xa = *colptr;
|
|
kusano |
7d535a |
|
|
kusano |
7d535a |
val = (double *) SUPERLU_MALLOC(*nonz * sizeof(double));
|
|
kusano |
7d535a |
row = (int *) SUPERLU_MALLOC(*nonz * sizeof(int));
|
|
kusano |
7d535a |
col = (int *) SUPERLU_MALLOC(*nonz * sizeof(int));
|
|
kusano |
7d535a |
|
|
kusano |
7d535a |
for (j = 0; j < *n; ++j) xa[j] = 0;
|
|
kusano |
7d535a |
|
|
kusano |
7d535a |
/* Read into the triplet array from a file */
|
|
kusano |
7d535a |
for (nnz = 0, nz = 0; nnz < *nonz; ++nnz) {
|
|
kusano |
7d535a |
scanf("%d%d%lf\n", &row[nz], &col[nz], &val[nz]);
|
|
kusano |
7d535a |
/* Change to 0-based indexing. */
|
|
kusano |
7d535a |
#if 0
|
|
kusano |
7d535a |
--row[nz];
|
|
kusano |
7d535a |
--col[nz];
|
|
kusano |
7d535a |
#endif
|
|
kusano |
7d535a |
if (row[nz] < 0 || row[nz] >= *m || col[nz] < 0 || col[nz] >= *n
|
|
kusano |
7d535a |
/*|| val[nz] == 0.*/) {
|
|
kusano |
7d535a |
fprintf(stderr, "nz %d, (%d, %d) = %e out of bound, removed\n",
|
|
kusano |
7d535a |
nz, row[nz], col[nz], val[nz]);
|
|
kusano |
7d535a |
exit(-1);
|
|
kusano |
7d535a |
} else {
|
|
kusano |
7d535a |
++xa[col[nz]];
|
|
kusano |
7d535a |
++nz;
|
|
kusano |
7d535a |
}
|
|
kusano |
7d535a |
}
|
|
kusano |
7d535a |
|
|
kusano |
7d535a |
*nonz = nz;
|
|
kusano |
7d535a |
|
|
kusano |
7d535a |
/* Initialize the array of column pointers */
|
|
kusano |
7d535a |
k = 0;
|
|
kusano |
7d535a |
jsize = xa[0];
|
|
kusano |
7d535a |
xa[0] = 0;
|
|
kusano |
7d535a |
for (j = 1; j < *n; ++j) {
|
|
kusano |
7d535a |
k += jsize;
|
|
kusano |
7d535a |
jsize = xa[j];
|
|
kusano |
7d535a |
xa[j] = k;
|
|
kusano |
7d535a |
}
|
|
kusano |
7d535a |
|
|
kusano |
7d535a |
/* Copy the triplets into the column oriented storage */
|
|
kusano |
7d535a |
for (nz = 0; nz < *nonz; ++nz) {
|
|
kusano |
7d535a |
j = col[nz];
|
|
kusano |
7d535a |
k = xa[j];
|
|
kusano |
7d535a |
asub[k] = row[nz];
|
|
kusano |
7d535a |
a[k] = val[nz];
|
|
kusano |
7d535a |
++xa[j];
|
|
kusano |
7d535a |
}
|
|
kusano |
7d535a |
|
|
kusano |
7d535a |
/* Reset the column pointers to the beginning of each column */
|
|
kusano |
7d535a |
for (j = *n; j > 0; --j)
|
|
kusano |
7d535a |
xa[j] = xa[j-1];
|
|
kusano |
7d535a |
xa[0] = 0;
|
|
kusano |
7d535a |
|
|
kusano |
7d535a |
SUPERLU_FREE(val);
|
|
kusano |
7d535a |
SUPERLU_FREE(row);
|
|
kusano |
7d535a |
SUPERLU_FREE(col);
|
|
kusano |
7d535a |
|
|
kusano |
7d535a |
#ifdef CHK_INPUT
|
|
kusano |
7d535a |
for (i = 0; i < *n; i++) {
|
|
kusano |
7d535a |
printf("Col %d, xa %d\n", i, xa[i]);
|
|
kusano |
7d535a |
for (k = xa[i]; k < xa[i+1]; k++)
|
|
kusano |
7d535a |
printf("%d\t%16.10f\n", asub[k], a[k]);
|
|
kusano |
7d535a |
}
|
|
kusano |
7d535a |
#endif
|
|
kusano |
7d535a |
|
|
kusano |
7d535a |
}
|
|
kusano |
7d535a |
|
|
kusano |
7d535a |
|
|
kusano |
7d535a |
void dreadrhs(int m, double *b)
|
|
kusano |
7d535a |
{
|
|
kusano |
7d535a |
FILE *fp, *fopen();
|
|
kusano |
7d535a |
int i, j;
|
|
kusano |
7d535a |
|
|
kusano |
7d535a |
if ( !(fp = fopen("b.dat", "r")) ) {
|
|
kusano |
7d535a |
fprintf(stderr, "dreadrhs: file does not exist\n");
|
|
kusano |
7d535a |
exit(-1);
|
|
kusano |
7d535a |
}
|
|
kusano |
7d535a |
for (i = 0; i < m; ++i)
|
|
kusano |
7d535a |
fscanf(fp, "%lf\n", &b[i]);
|
|
kusano |
7d535a |
/*fscanf(fp, "%d%lf\n", &j, &b[i]);*/
|
|
kusano |
7d535a |
/* readpair_(j, &b[i]);*/
|
|
kusano |
7d535a |
fclose(fp);
|
|
kusano |
7d535a |
}
|
|
kusano |
7d535a |
|
|
kusano |
7d535a |
|