1 /*! \file
2 Copyright (c) 2003, The Regents of the University of California, through
3 Lawrence Berkeley National Laboratory (subject to receipt of any required
4 approvals from U.S. Dept. of Energy)
5
6 All rights reserved.
7
8 The source code is distributed under BSD license, see the file License.txt
9 at the top-level directory.
10 */
11
12 /*! @file dreadtriple.c
13 * \brief Read a matrix stored in triplet (coordinate) format
14 *
15 * <pre>
16 * -- SuperLU routine (version 4.0) --
17 * Lawrence Berkeley National Laboratory.
18 * June 30, 2009
19 * </pre>
20 */
21
22 #include "slu_ddefs.h"
23
24
25 void
dreadtriple(int * m,int * n,int * nonz,double ** nzval,int ** rowind,int ** colptr)26 dreadtriple(int *m, int *n, int *nonz,
27 double **nzval, int **rowind, int **colptr)
28 {
29 /*
30 * Output parameters
31 * =================
32 * (a,asub,xa): asub[*] contains the row subscripts of nonzeros
33 * in columns of matrix A; a[*] the numerical values;
34 * row i of A is given by a[k],k=xa[i],...,xa[i+1]-1.
35 *
36 */
37 int j, k, jsize, nnz, nz;
38 double *a, *val;
39 int *asub, *xa, *row, *col;
40 int zero_base = 0;
41
42 /* Matrix format:
43 * First line: #rows, #cols, #non-zero
44 * Triplet in the rest of lines:
45 * row, col, value
46 */
47
48 scanf("%d%d", n, nonz);
49 *m = *n;
50 printf("m %d, n %d, nonz %d\n", *m, *n, *nonz);
51 dallocateA(*n, *nonz, nzval, rowind, colptr); /* Allocate storage */
52 a = *nzval;
53 asub = *rowind;
54 xa = *colptr;
55
56 val = (double *) SUPERLU_MALLOC(*nonz * sizeof(double));
57 row = (int *) SUPERLU_MALLOC(*nonz * sizeof(int));
58 col = (int *) SUPERLU_MALLOC(*nonz * sizeof(int));
59
60 for (j = 0; j < *n; ++j) xa[j] = 0;
61
62 /* Read into the triplet array from a file */
63 for (nnz = 0, nz = 0; nnz < *nonz; ++nnz) {
64 scanf("%d%d%lf\n", &row[nz], &col[nz], &val[nz]);
65
66 if ( nnz == 0 ) { /* first nonzero */
67 if ( row[0] == 0 || col[0] == 0 ) {
68 zero_base = 1;
69 printf("triplet file: row/col indices are zero-based.\n");
70 } else
71 printf("triplet file: row/col indices are one-based.\n");
72 }
73
74 if ( !zero_base ) {
75 /* Change to 0-based indexing. */
76 --row[nz];
77 --col[nz];
78 }
79
80 if (row[nz] < 0 || row[nz] >= *m || col[nz] < 0 || col[nz] >= *n
81 /*|| val[nz] == 0.*/) {
82 fprintf(stderr, "nz %d, (%d, %d) = %e out of bound, removed\n",
83 nz, row[nz], col[nz], val[nz]);
84 exit(-1);
85 } else {
86 ++xa[col[nz]];
87 ++nz;
88 }
89 }
90
91 *nonz = nz;
92
93 /* Initialize the array of column pointers */
94 k = 0;
95 jsize = xa[0];
96 xa[0] = 0;
97 for (j = 1; j < *n; ++j) {
98 k += jsize;
99 jsize = xa[j];
100 xa[j] = k;
101 }
102
103 /* Copy the triplets into the column oriented storage */
104 for (nz = 0; nz < *nonz; ++nz) {
105 j = col[nz];
106 k = xa[j];
107 asub[k] = row[nz];
108 a[k] = val[nz];
109 ++xa[j];
110 }
111
112 /* Reset the column pointers to the beginning of each column */
113 for (j = *n; j > 0; --j)
114 xa[j] = xa[j-1];
115 xa[0] = 0;
116
117 SUPERLU_FREE(val);
118 SUPERLU_FREE(row);
119 SUPERLU_FREE(col);
120
121 #ifdef CHK_INPUT
122 {
123 int i;
124 for (i = 0; i < *n; i++) {
125 printf("Col %d, xa %d\n", i, xa[i]);
126 for (k = xa[i]; k < xa[i+1]; k++)
127 printf("%d\t%16.10f\n", asub[k], a[k]);
128 }
129 }
130 #endif
131
132 }
133
134
dreadrhs(int m,double * b)135 void dreadrhs(int m, double *b)
136 {
137 FILE *fp, *fopen();
138 int i;
139 /*int j;*/
140
141 if ( !(fp = fopen("b.dat", "r")) ) {
142 fprintf(stderr, "dreadrhs: file does not exist\n");
143 exit(-1);
144 }
145 for (i = 0; i < m; ++i)
146 fscanf(fp, "%lf\n", &b[i]);
147
148 /* readpair_(j, &b[i]);*/
149 fclose(fp);
150 }
151