1/* $Id$Revision: */
2/* vim:set shiftwidth=4 ts=8: */
3
4/*************************************************************************
5 * Copyright (c) 2011 AT&T Intellectual Property
6 * All rights reserved. This program and the accompanying materials
7 * are made available under the terms of the Eclipse Public License v1.0
8 * which accompanies this distribution, and is available at
9 * http://www.eclipse.org/legal/epl-v10.html
10 *
11 * Contributors: See CVS logs. Details at http://www.graphviz.org/
12 *************************************************************************/
13
14#include "SparseMatrix.h"
15#include "mmio.h"
16#include "matrix_market.h"
17#include "memory.h"
18#include "assert.h"
19#define MALLOC gmalloc
20#define REALLOC grealloc
21#define FREE free
22#define MEMCPY memcpy
23
24int mm_get_type(MM_typecode typecode)
25{
26 if (mm_is_complex(typecode)) {
27 return MATRIX_TYPE_COMPLEX;
28 } else if (mm_is_real(typecode)) {
29 return MATRIX_TYPE_REAL;
30 } else if (mm_is_integer(typecode)) {
31 return MATRIX_TYPE_INTEGER;
32 } else if (mm_is_pattern(typecode)) {
33 return MATRIX_TYPE_PATTERN;
34 }
35 return MATRIX_TYPE_UNKNOWN;
36}
37
38static void set_mm_typecode(int type, MM_typecode * typecode)
39{
40 switch (type) {
41 case MATRIX_TYPE_COMPLEX:
42 mm_set_complex(typecode);
43 break;
44 case MATRIX_TYPE_REAL:
45 mm_set_real(typecode);
46 break;
47 case MATRIX_TYPE_INTEGER:
48 mm_set_integer(typecode);
49 break;
50 case MATRIX_TYPE_PATTERN:
51 mm_set_pattern(typecode);
52 break;
53 default:
54 break;
55 }
56}
57
58
59
60
61SparseMatrix SparseMatrix_import_matrix_market(FILE * f, int format)
62{
63 int ret_code, type;
64 MM_typecode matcode;
65 real *val = NULL, *v;
66 int *vali = NULL, i, m, n, *I = NULL, *J = NULL, nz;
67 void *vp = NULL;
68 SparseMatrix A = NULL;
69 int nzold;
70 int c;
71
72 if ((c = fgetc(f)) != '%') {
73 ungetc(c, f);
74 return NULL;
75 }
76 ungetc(c, f);
77 if (mm_read_banner(f, &matcode) != 0) {
78#ifdef DEBUG
79 printf("Could not process Matrix Market banner.\n");
80#endif
81 return NULL;
82 }
83
84
85 /* This is how one can screen matrix types if their application */
86 /* only supports a subset of the Matrix Market data types. */
87
88 if (!mm_is_matrix(matcode) || !mm_is_sparse(matcode)) {
89 assert(0);
90 /*
91 printf("Sorry, this application does not support dense matrix");
92 printf("Market Market type: [%s]\n", mm_typecode_to_str(matcode));
93 */
94 return NULL;
95 }
96
97 /* find out size of sparse matrix .... */
98 if ((ret_code = mm_read_mtx_crd_size(f, &m, &n, &nz)) != 0) {
99 assert(0);
100 return NULL;
101 }
102 /* reseve memory for matrices */
103
104 I = MALLOC(nz * sizeof(int));
105 J = MALLOC(nz * sizeof(int));
106
107
108
109 switch (format) {
110 case FORMAT_CSC:
111 assert(0); /* not supported yet */
112 break;
113 case FORMAT_CSR:
114 case FORMAT_COORD:
115
116 /* NOTE: when reading in doubles, ANSI C requires the use of the "l" */
117 /* specifier as in "%lg", "%lf", "%le", otherwise errors will occur */
118 /* (ANSI C X3.159-1989, Sec. 4.9.6.2, p. 136 lines 13-15) */
119 type = mm_get_type(matcode);
120 switch (type) {
121 case MATRIX_TYPE_REAL:
122 val = (real *) malloc(nz * sizeof(real));
123 for (i = 0; i < nz; i++) {
124 fscanf(f, "%d %d %lg\n", &I[i], &J[i], &val[i]);
125 I[i]--; /* adjust from 1-based to 0-based */
126 J[i]--;
127 }
128 if (mm_is_symmetric(matcode)) {
129 I = REALLOC(I, 2 * sizeof(int) * nz);
130 J = REALLOC(J, 2 * sizeof(int) * nz);
131 val = REALLOC(val, 2 * sizeof(real) * nz);
132 nzold = nz;
133 for (i = 0; i < nzold; i++) {
134 if (I[i] != J[i]) {
135 I[nz] = J[i];
136 J[nz] = I[i];
137 val[nz++] = val[i];
138 }
139 }
140 } else if (mm_is_skew(matcode)) {
141 I = REALLOC(I, 2 * sizeof(int) * nz);
142 J = REALLOC(J, 2 * sizeof(int) * nz);
143 val = REALLOC(val, 2 * sizeof(real) * nz);
144 vp = (void *) val;
145 nzold = nz;
146 for (i = 0; i < nzold; i++) {
147 assert(I[i] != J[i]); /* skew symm has no diag */
148 I[nz] = J[i];
149 J[nz] = I[i];
150 val[nz++] = -val[i];
151 }
152 } else {
153 assert(!mm_is_hermitian(matcode));
154 }
155 vp = (void *) val;
156 break;
157 case MATRIX_TYPE_INTEGER:
158 vali = (int *) malloc(nz * sizeof(int));
159 for (i = 0; i < nz; i++) {
160 fscanf(f, "%d %d %d\n", &I[i], &J[i], &vali[i]);
161 I[i]--; /* adjust from 1-based to 0-based */
162 J[i]--;
163 }
164 if (mm_is_symmetric(matcode)) {
165 I = REALLOC(I, 2 * sizeof(int) * nz);
166 J = REALLOC(J, 2 * sizeof(int) * nz);
167 vali = REALLOC(vali, 2 * sizeof(int) * nz);
168 nzold = nz;
169 for (i = 0; i < nzold; i++) {
170 if (I[i] != J[i]) {
171 I[nz] = J[i];
172 J[nz] = I[i];
173 vali[nz++] = vali[i];
174 }
175 }
176 } else if (mm_is_skew(matcode)) {
177 I = REALLOC(I, 2 * sizeof(int) * nz);
178 J = REALLOC(J, 2 * sizeof(int) * nz);
179 vali = REALLOC(vali, 2 * sizeof(int) * nz);
180 vp = (void *) val;
181 nzold = nz;
182 for (i = 0; i < nzold; i++) {
183 assert(I[i] != J[i]); /* skew symm has no diag */
184 I[nz] = J[i];
185 J[nz] = I[i];
186 vali[nz++] = -vali[i];
187 }
188 } else {
189 assert(!mm_is_hermitian(matcode));
190 }
191 vp = (void *) vali;
192 break;
193 case MATRIX_TYPE_PATTERN:
194 for (i = 0; i < nz; i++) {
195 fscanf(f, "%d %d\n", &I[i], &J[i]);
196 I[i]--; /* adjust from 1-based to 0-based */
197 J[i]--;
198 }
199 if (mm_is_symmetric(matcode) || mm_is_skew(matcode)) {
200 I = REALLOC(I, 2 * sizeof(int) * nz);
201 J = REALLOC(J, 2 * sizeof(int) * nz);
202 nzold = nz;
203 for (i = 0; i < nzold; i++) {
204 if (I[i] != J[i]) {
205 I[nz] = J[i];
206 J[nz++] = I[i];
207 }
208 }
209 } else {
210 assert(!mm_is_hermitian(matcode));
211 }
212 break;
213 case MATRIX_TYPE_COMPLEX:
214 val = (real *) malloc(2 * nz * sizeof(real));
215 v = val;
216 for (i = 0; i < nz; i++) {
217 fscanf(f, "%d %d %lg %lg\n", &I[i], &J[i], &v[0], &v[1]);
218 v += 2;
219 I[i]--; /* adjust from 1-based to 0-based */
220 J[i]--;
221 }
222 if (mm_is_symmetric(matcode)) {
223 I = REALLOC(I, 2 * sizeof(int) * nz);
224 J = REALLOC(J, 2 * sizeof(int) * nz);
225 val = REALLOC(val, 4 * sizeof(real) * nz);
226 nzold = nz;
227 for (i = 0; i < nzold; i++) {
228 if (I[i] != J[i]) {
229 I[nz] = J[i];
230 J[nz] = I[i];
231 val[2 * nz] = val[2 * i];
232 val[2 * nz + 1] = val[2 * i + 1];
233 nz++;
234 }
235 }
236 } else if (mm_is_skew(matcode)) {
237 I = REALLOC(I, 2 * sizeof(int) * nz);
238 J = REALLOC(J, 2 * sizeof(int) * nz);
239 val = REALLOC(val, 4 * sizeof(real) * nz);
240 vp = (void *) val;
241 nzold = nz;
242 for (i = 0; i < nzold; i++) {
243 assert(I[i] != J[i]); /* skew symm has no diag */
244 I[nz] = J[i];
245 J[nz] = I[i];
246 val[2 * nz] = -val[2 * i];
247 val[2 * nz + 1] = -val[2 * i + 1];
248 nz++;
249
250 }
251 } else if (mm_is_hermitian(matcode)) {
252 I = REALLOC(I, 2 * sizeof(int) * nz);
253 J = REALLOC(J, 2 * sizeof(int) * nz);
254 val = REALLOC(val, 4 * sizeof(real) * nz);
255 vp = (void *) val;
256 nzold = nz;
257 for (i = 0; i < nzold; i++) {
258 if (I[i] != J[i]) {
259 I[nz] = J[i];
260 J[nz] = I[i];
261 val[2 * nz] = val[2 * i];
262 val[2 * nz + 1] = -val[2 * i + 1];
263 nz++;
264 }
265 }
266 }
267 vp = (void *) val;
268 break;
269 default:
270 return 0;
271 }
272
273 if (format == FORMAT_CSR) {
274 A = SparseMatrix_from_coordinate_arrays(nz, m, n, I, J, vp,
275 type, sizeof(real));
276 } else {
277 A = SparseMatrix_new(m, n, 1, type, FORMAT_COORD);
278 A = SparseMatrix_coordinate_form_add_entries(A, nz, I, J, vp);
279 }
280 break;
281 default:
282 A = NULL;
283 }
284 FREE(I);
285 FREE(J);
286 FREE(val);
287
288 if (mm_is_symmetric(matcode)) {
289 SparseMatrix_set_symmetric(A);
290 SparseMatrix_set_pattern_symmetric(A);
291 } else if (mm_is_skew(matcode)) {
292 SparseMatrix_set_skew(A);
293 } else if (mm_is_hermitian(matcode)) {
294 SparseMatrix_set_hemitian(A);
295 }
296
297
298 return A;
299}
300
301
302static void mm_write_comment(FILE * file, char *comment)
303{
304 char percent[2] = "%";
305 fprintf(file, "%s %s\n", percent, comment);
306}
307
308void SparseMatrix_export_matrix_market(FILE * file, SparseMatrix A,
309 char *comment)
310{
311 MM_typecode matcode;
312
313 mm_initialize_typecode(&matcode);
314 mm_set_matrix(&matcode);
315 mm_set_sparse(&matcode);
316 mm_set_general(&matcode);
317 set_mm_typecode(A->type, &matcode);
318
319 mm_write_banner(file, matcode);
320 mm_write_comment(file, comment);
321
322 SparseMatrix_export(file, A);
323
324}
325