1 |
|
2 |
/******************************************************* |
3 |
* |
4 |
* Copyright (c) 2003-2008 by University of Queensland |
5 |
* Earth Systems Science Computational Center (ESSCC) |
6 |
* http://www.uq.edu.au/esscc |
7 |
* |
8 |
* Primary Business: Queensland, Australia |
9 |
* Licensed under the Open Software License version 3.0 |
10 |
* http://www.opensource.org/licenses/osl-3.0.php |
11 |
* |
12 |
*******************************************************/ |
13 |
|
14 |
|
15 |
/* |
16 |
* Matrix Market I/O library for ANSI C |
17 |
* |
18 |
* See http://math.nist.gov/MatrixMarket for details. |
19 |
* |
20 |
* (Version 1.01, 5/2003) |
21 |
*/ |
22 |
|
23 |
#include "Common.h" |
24 |
#include "mmio.h" |
25 |
|
26 |
#define FSCANF_CHECK(scan_ret, reason) { if (scan_ret == EOF) perror(reason); return -1; } |
27 |
|
28 |
#include <string.h> |
29 |
#include <ctype.h> |
30 |
|
31 |
|
32 |
|
33 |
int mm_read_unsymmetric_sparse(const char *fname, int *M_, int *N_, int *nz_, |
34 |
double **val_, int **I_, int **J_) |
35 |
{ |
36 |
FILE *f; |
37 |
MM_typecode matcode; |
38 |
int M, N, nz; |
39 |
int i; |
40 |
double *val; |
41 |
int *Ip, *Jp; |
42 |
|
43 |
if ((f = fopen(fname, "r")) == NULL) |
44 |
return -1; |
45 |
|
46 |
|
47 |
if (mm_read_banner(f, &matcode) != 0) |
48 |
{ |
49 |
printf("mm_read_unsymetric: Could not process Matrix Market banner "); |
50 |
printf(" in file [%s]\n", fname); |
51 |
return -1; |
52 |
} |
53 |
|
54 |
|
55 |
|
56 |
if ( !(mm_is_real(matcode) && mm_is_matrix(matcode) && |
57 |
mm_is_sparse(matcode))) |
58 |
{ |
59 |
fprintf(stderr, "Sorry, this application does not support "); |
60 |
fprintf(stderr, "Market Market type: [%s]\n", |
61 |
mm_typecode_to_str(matcode)); |
62 |
return -1; |
63 |
} |
64 |
|
65 |
/* find out size of sparse matrix: M, N, nz .... */ |
66 |
|
67 |
if (mm_read_mtx_crd_size(f, &M, &N, &nz) !=0) |
68 |
{ |
69 |
fprintf(stderr, "read_unsymmetric_sparse(): could not parse matrix size.\n"); |
70 |
return -1; |
71 |
} |
72 |
|
73 |
*M_ = M; |
74 |
*N_ = N; |
75 |
*nz_ = nz; |
76 |
|
77 |
/* reseve memory for matrices */ |
78 |
|
79 |
Ip = MEMALLOC(nz, int); |
80 |
Jp = MEMALLOC(nz, int); |
81 |
val = MEMALLOC(nz, double); |
82 |
|
83 |
*val_ = val; |
84 |
*I_ = Ip; |
85 |
*J_ = Jp; |
86 |
|
87 |
/* NOTE: when reading in doubles, ANSI C requires the use of the "l" */ |
88 |
/* specifier as in "%lg", "%lf", "%le", otherwise errors will occur */ |
89 |
/* (ANSI C X3.159-1989, Sec. 4.9.6.2, p. 136 lines 13-15) */ |
90 |
|
91 |
for (i=0; i<nz; i++) |
92 |
{ |
93 |
int scan_ret = fscanf(f, "%d %d %lg\n", &Ip[i], &Jp[i], &val[i]); |
94 |
FSCANF_CHECK(scan_ret, "fscanf: mm_read_unsymmetric_sparse"); |
95 |
Ip[i]--; /* adjust from 1-based to 0-based */ |
96 |
Jp[i]--; |
97 |
} |
98 |
fclose(f); |
99 |
|
100 |
return 0; |
101 |
} |
102 |
|
103 |
int mm_is_valid(MM_typecode matcode) |
104 |
{ |
105 |
if (!mm_is_matrix(matcode)) return 0; |
106 |
if (mm_is_dense(matcode) && mm_is_pattern(matcode)) return 0; |
107 |
if (mm_is_real(matcode) && mm_is_hermitian(matcode)) return 0; |
108 |
if (mm_is_pattern(matcode) && (mm_is_hermitian(matcode) || |
109 |
mm_is_skew(matcode))) return 0; |
110 |
return 1; |
111 |
} |
112 |
|
113 |
int mm_read_banner(FILE *f, MM_typecode *matcode) |
114 |
{ |
115 |
char line[MM_MAX_LINE_LENGTH]; |
116 |
char banner[MM_MAX_TOKEN_LENGTH]; |
117 |
char mtx[MM_MAX_TOKEN_LENGTH]; |
118 |
char crd[MM_MAX_TOKEN_LENGTH]; |
119 |
char data_type[MM_MAX_TOKEN_LENGTH]; |
120 |
char storage_scheme[MM_MAX_TOKEN_LENGTH]; |
121 |
char *p; |
122 |
|
123 |
|
124 |
mm_clear_typecode(matcode); |
125 |
|
126 |
if (fgets(line, MM_MAX_LINE_LENGTH, f) == NULL) |
127 |
return MM_PREMATURE_EOF; |
128 |
|
129 |
if (sscanf(line, "%s %s %s %s %s", banner, mtx, crd, data_type, |
130 |
storage_scheme) != 5) |
131 |
return MM_PREMATURE_EOF; |
132 |
|
133 |
/* convert to lower case */ |
134 |
for (p=mtx; *p!='\0'; *p= (char) tolower(*p),p++); |
135 |
for (p=crd; *p!='\0'; *p= (char) tolower(*p),p++); |
136 |
for (p=data_type; *p!='\0'; *p= (char) tolower(*p),p++); |
137 |
for (p=storage_scheme; *p!='\0'; *p= (char) tolower(*p),p++); |
138 |
|
139 |
/* check for banner */ |
140 |
if (strncmp(banner, MatrixMarketBanner, strlen(MatrixMarketBanner)) != 0) |
141 |
return MM_NO_HEADER; |
142 |
|
143 |
/* first field should be "mtx" */ |
144 |
if (strcmp(mtx, MM_MTX_STR) != 0) |
145 |
return MM_UNSUPPORTED_TYPE; |
146 |
mm_set_matrix(matcode); |
147 |
|
148 |
|
149 |
/* second field describes whether this is a sparse matrix (in coordinate |
150 |
storgae) or a dense array */ |
151 |
|
152 |
|
153 |
if (strcmp(crd, MM_SPARSE_STR) == 0) |
154 |
mm_set_sparse(matcode); |
155 |
else |
156 |
if (strcmp(crd, MM_DENSE_STR) == 0) |
157 |
mm_set_dense(matcode); |
158 |
else |
159 |
return MM_UNSUPPORTED_TYPE; |
160 |
|
161 |
|
162 |
/* third field */ |
163 |
|
164 |
if (strcmp(data_type, MM_REAL_STR) == 0) |
165 |
mm_set_real(matcode); |
166 |
else |
167 |
if (strcmp(data_type, MM_COMPLEX_STR) == 0) |
168 |
mm_set_complex(matcode); |
169 |
else |
170 |
if (strcmp(data_type, MM_PATTERN_STR) == 0) |
171 |
mm_set_pattern(matcode); |
172 |
else |
173 |
if (strcmp(data_type, MM_INT_STR) == 0) |
174 |
mm_set_integer(matcode); |
175 |
else |
176 |
return MM_UNSUPPORTED_TYPE; |
177 |
|
178 |
|
179 |
/* fourth field */ |
180 |
|
181 |
if (strcmp(storage_scheme, MM_GENERAL_STR) == 0) |
182 |
mm_set_general(matcode); |
183 |
else |
184 |
if (strcmp(storage_scheme, MM_SYMM_STR) == 0) |
185 |
mm_set_symmetric(matcode); |
186 |
else |
187 |
if (strcmp(storage_scheme, MM_HERM_STR) == 0) |
188 |
mm_set_hermitian(matcode); |
189 |
else |
190 |
if (strcmp(storage_scheme, MM_SKEW_STR) == 0) |
191 |
mm_set_skew(matcode); |
192 |
else |
193 |
return MM_UNSUPPORTED_TYPE; |
194 |
|
195 |
|
196 |
return 0; |
197 |
} |
198 |
|
199 |
int mm_write_mtx_crd_size(FILE *f, int M, int N, int nz) |
200 |
{ |
201 |
if (fprintf(f, "%d %d %d\n", M, N, nz) < 0) |
202 |
return MM_COULD_NOT_WRITE_FILE; |
203 |
else |
204 |
return 0; |
205 |
} |
206 |
|
207 |
int mm_read_mtx_crd_size(FILE *f, int *M, int *N, int *nz ) |
208 |
{ |
209 |
char line[MM_MAX_LINE_LENGTH]; |
210 |
int num_items_read; |
211 |
|
212 |
/* set return null parameter values, in case we exit with errors */ |
213 |
*M = *N = *nz = 0; |
214 |
|
215 |
/* now continue scanning until you reach the end-of-comments */ |
216 |
do |
217 |
{ |
218 |
if (fgets(line,MM_MAX_LINE_LENGTH,f) == NULL) |
219 |
return MM_PREMATURE_EOF; |
220 |
}while (line[0] == '%'); |
221 |
|
222 |
/* line[] is either blank or has M,N, nz */ |
223 |
if (sscanf(line, "%d %d %d", M, N, nz) == 3) |
224 |
return 0; |
225 |
|
226 |
else |
227 |
do |
228 |
{ |
229 |
num_items_read = fscanf(f, "%d %d %d", M, N, nz); |
230 |
if (num_items_read == EOF) return MM_PREMATURE_EOF; |
231 |
} |
232 |
while (num_items_read != 3); |
233 |
|
234 |
return 0; |
235 |
} |
236 |
|
237 |
|
238 |
int mm_read_mtx_array_size(FILE *f, int *M, int *N) |
239 |
{ |
240 |
char line[MM_MAX_LINE_LENGTH]; |
241 |
int num_items_read; |
242 |
/* set return null parameter values, in case we exit with errors */ |
243 |
*M = *N = 0; |
244 |
|
245 |
/* now continue scanning until you reach the end-of-comments */ |
246 |
do |
247 |
{ |
248 |
if (fgets(line,MM_MAX_LINE_LENGTH,f) == NULL) |
249 |
return MM_PREMATURE_EOF; |
250 |
}while (line[0] == '%'); |
251 |
|
252 |
/* line[] is either blank or has M,N, nz */ |
253 |
if (sscanf(line, "%d %d", M, N) == 2) |
254 |
return 0; |
255 |
|
256 |
else /* we have a blank line */ |
257 |
do |
258 |
{ |
259 |
num_items_read = fscanf(f, "%d %d", M, N); |
260 |
if (num_items_read == EOF) return MM_PREMATURE_EOF; |
261 |
} |
262 |
while (num_items_read != 2); |
263 |
|
264 |
return 0; |
265 |
} |
266 |
|
267 |
int mm_write_mtx_array_size(FILE *f, int M, int N) |
268 |
{ |
269 |
if (fprintf(f, "%d %d\n", M, N) < 0) |
270 |
return MM_COULD_NOT_WRITE_FILE; |
271 |
else |
272 |
return 0; |
273 |
} |
274 |
|
275 |
|
276 |
|
277 |
/*-------------------------------------------------------------------------*/ |
278 |
|
279 |
/******************************************************************/ |
280 |
/* use when Ip[], Jp[], and val[] are already allocated */ |
281 |
/******************************************************************/ |
282 |
|
283 |
int mm_read_mtx_crd_data(FILE *f, int M, int N, int nz, int Ip[], int Jp[], |
284 |
double val[], MM_typecode matcode) |
285 |
{ |
286 |
int i; |
287 |
if (mm_is_complex(matcode)) |
288 |
{ |
289 |
for (i=0; i<nz; i++) |
290 |
if (fscanf(f, "%d %d %lg %lg", &Ip[i], &Jp[i], &val[2*i], &val[2*i+1]) |
291 |
!= 4) return MM_PREMATURE_EOF; |
292 |
} |
293 |
else if (mm_is_real(matcode)) |
294 |
{ |
295 |
for (i=0; i<nz; i++) |
296 |
{ |
297 |
if (fscanf(f, "%d %d %lg\n", &Ip[i], &Jp[i], &val[i]) |
298 |
!= 3) return MM_PREMATURE_EOF; |
299 |
|
300 |
} |
301 |
} |
302 |
|
303 |
else if (mm_is_pattern(matcode)) |
304 |
{ |
305 |
for (i=0; i<nz; i++) |
306 |
if (fscanf(f, "%d %d", &Ip[i], &Jp[i]) |
307 |
!= 2) return MM_PREMATURE_EOF; |
308 |
} |
309 |
else |
310 |
return MM_UNSUPPORTED_TYPE; |
311 |
|
312 |
return 0; |
313 |
|
314 |
} |
315 |
|
316 |
int mm_read_mtx_crd_entry(FILE *f, int *Ip, int *Jp, |
317 |
double *real, double *imag, MM_typecode matcode) |
318 |
{ |
319 |
if (mm_is_complex(matcode)) |
320 |
{ |
321 |
if (fscanf(f, "%d %d %lg %lg", Ip, Jp, real, imag) |
322 |
!= 4) return MM_PREMATURE_EOF; |
323 |
} |
324 |
else if (mm_is_real(matcode)) |
325 |
{ |
326 |
if (fscanf(f, "%d %d %lg\n", Ip, Jp, real) |
327 |
!= 3) return MM_PREMATURE_EOF; |
328 |
|
329 |
} |
330 |
|
331 |
else if (mm_is_pattern(matcode)) |
332 |
{ |
333 |
if (fscanf(f, "%d %d", Ip, Jp) != 2) return MM_PREMATURE_EOF; |
334 |
} |
335 |
else |
336 |
return MM_UNSUPPORTED_TYPE; |
337 |
|
338 |
return 0; |
339 |
|
340 |
} |
341 |
|
342 |
|
343 |
/************************************************************************ |
344 |
mm_read_mtx_crd() fills M, N, nz, array of values, and return |
345 |
type code, e.g. 'MCRS' |
346 |
|
347 |
if matrix is complex, values[] is of size 2*nz, |
348 |
(nz pairs of real/imaginary values) |
349 |
************************************************************************/ |
350 |
|
351 |
int mm_read_mtx_crd(char *fname, int *M, int *N, int *nz, int **Ip, int **Jp, |
352 |
double **val, MM_typecode *matcode) |
353 |
{ |
354 |
int ret_code; |
355 |
FILE *f; |
356 |
|
357 |
if (strcmp(fname, "stdin") == 0) f=stdin; |
358 |
else |
359 |
if ((f = fopen(fname, "r")) == NULL) |
360 |
return MM_COULD_NOT_READ_FILE; |
361 |
|
362 |
|
363 |
if ((ret_code = mm_read_banner(f, matcode)) != 0) |
364 |
return ret_code; |
365 |
|
366 |
if (!(mm_is_valid(*matcode) && mm_is_sparse(*matcode) && |
367 |
mm_is_matrix(*matcode))) |
368 |
return MM_UNSUPPORTED_TYPE; |
369 |
|
370 |
if ((ret_code = mm_read_mtx_crd_size(f, M, N, nz)) != 0) |
371 |
return ret_code; |
372 |
|
373 |
|
374 |
*Ip = MEMALLOC(*nz, int); |
375 |
*Jp = MEMALLOC(*nz, int); |
376 |
*val = NULL; |
377 |
|
378 |
if (mm_is_complex(*matcode)) |
379 |
{ |
380 |
*val = MEMALLOC(*nz * 2, double); |
381 |
ret_code = mm_read_mtx_crd_data(f, *M, *N, *nz, *Ip, *Jp, *val, |
382 |
*matcode); |
383 |
if (ret_code != 0) return ret_code; |
384 |
} |
385 |
else if (mm_is_real(*matcode)) |
386 |
{ |
387 |
*val = MEMALLOC(*nz, double); |
388 |
ret_code = mm_read_mtx_crd_data(f, *M, *N, *nz, *Ip, *Jp, *val, |
389 |
*matcode); |
390 |
if (ret_code != 0) return ret_code; |
391 |
} |
392 |
|
393 |
else if (mm_is_pattern(*matcode)) |
394 |
{ |
395 |
ret_code = mm_read_mtx_crd_data(f, *M, *N, *nz, *Ip, *Jp, *val, |
396 |
*matcode); |
397 |
if (ret_code != 0) return ret_code; |
398 |
} |
399 |
|
400 |
if (f != stdin) fclose(f); |
401 |
return 0; |
402 |
} |
403 |
|
404 |
int mm_write_banner(FILE *f, MM_typecode matcode) |
405 |
{ |
406 |
int ret_code = fprintf(f, "%s %s\n", MatrixMarketBanner, mm_typecode_to_str(matcode)); |
407 |
if (ret_code < 0 ) |
408 |
return MM_COULD_NOT_WRITE_FILE; |
409 |
else |
410 |
return 0; |
411 |
} |
412 |
|
413 |
int mm_write_mtx_crd(char fname[], int M, int N, int nz, int Ip[], int Jp[], |
414 |
double val[], MM_typecode matcode) |
415 |
{ |
416 |
FILE *f; |
417 |
int i; |
418 |
|
419 |
if (strcmp(fname, "stdout") == 0) |
420 |
f = stdout; |
421 |
else |
422 |
if ((f = fopen(fname, "w")) == NULL) |
423 |
return MM_COULD_NOT_WRITE_FILE; |
424 |
|
425 |
/* print banner followed by typecode */ |
426 |
fprintf(f, "%s ", MatrixMarketBanner); |
427 |
fprintf(f, "%s\n", mm_typecode_to_str(matcode)); |
428 |
|
429 |
/* print matrix sizes and nonzeros */ |
430 |
fprintf(f, "%d %d %d\n", M, N, nz); |
431 |
|
432 |
/* print values */ |
433 |
if (mm_is_pattern(matcode)) |
434 |
for (i=0; i<nz; i++) |
435 |
fprintf(f, "%d %d\n", Ip[i], Jp[i]); |
436 |
else |
437 |
if (mm_is_real(matcode)) |
438 |
for (i=0; i<nz; i++) |
439 |
fprintf(f, "%d %d %20.16g\n", Ip[i], Jp[i], val[i]); |
440 |
else |
441 |
if (mm_is_complex(matcode)) |
442 |
for (i=0; i<nz; i++) |
443 |
fprintf(f, "%d %d %20.16g %20.16g\n", Ip[i], Jp[i], val[2*i], |
444 |
val[2*i+1]); |
445 |
else |
446 |
{ |
447 |
if (f != stdout) fclose(f); |
448 |
return MM_UNSUPPORTED_TYPE; |
449 |
} |
450 |
|
451 |
if (f !=stdout) fclose(f); |
452 |
|
453 |
return 0; |
454 |
} |
455 |
|
456 |
|
457 |
char *mm_typecode_to_str(MM_typecode matcode) |
458 |
{ |
459 |
static char buffer[MM_MAX_LINE_LENGTH]; |
460 |
char *types[4]; |
461 |
|
462 |
/* check for MTX type */ |
463 |
if (mm_is_matrix(matcode)) |
464 |
types[0] = MM_MTX_STR; |
465 |
else |
466 |
return NULL; |
467 |
|
468 |
/* check for CRD or ARR matrix */ |
469 |
if (mm_is_sparse(matcode)) |
470 |
types[1] = MM_SPARSE_STR; |
471 |
else |
472 |
if (mm_is_dense(matcode)) |
473 |
types[1] = MM_DENSE_STR; |
474 |
else |
475 |
return NULL; |
476 |
|
477 |
/* check for element data type */ |
478 |
if (mm_is_real(matcode)) |
479 |
types[2] = MM_REAL_STR; |
480 |
else |
481 |
if (mm_is_complex(matcode)) |
482 |
types[2] = MM_COMPLEX_STR; |
483 |
else |
484 |
if (mm_is_pattern(matcode)) |
485 |
types[2] = MM_PATTERN_STR; |
486 |
else |
487 |
if (mm_is_integer(matcode)) |
488 |
types[2] = MM_INT_STR; |
489 |
else |
490 |
return NULL; |
491 |
|
492 |
|
493 |
/* check for symmetry type */ |
494 |
if (mm_is_general(matcode)) |
495 |
types[3] = MM_GENERAL_STR; |
496 |
else |
497 |
if (mm_is_symmetric(matcode)) |
498 |
types[3] = MM_SYMM_STR; |
499 |
else |
500 |
if (mm_is_hermitian(matcode)) |
501 |
types[3] = MM_HERM_STR; |
502 |
else |
503 |
if (mm_is_skew(matcode)) |
504 |
types[3] = MM_SKEW_STR; |
505 |
else |
506 |
return NULL; |
507 |
|
508 |
sprintf(buffer,"%s %s %s %s", types[0], types[1], types[2], types[3]); |
509 |
return & buffer[0]; |
510 |
|
511 |
} |
512 |
|
513 |
/* |
514 |
* $Log$ |
515 |
* Revision 1.1 2004/10/26 06:53:59 jgs |
516 |
* Initial revision |
517 |
* |
518 |
* Revision 1.1 2004/07/02 00:48:35 gross |
519 |
* matrix market io function added |
520 |
* |
521 |
* |
522 |
*/ |