/[escript]/trunk/paso/src/mmio.c
ViewVC logotype

Contents of /trunk/paso/src/mmio.c

Parent Directory Parent Directory | Revision Log Revision Log


Revision 2548 - (show annotations)
Mon Jul 20 06:20:06 2009 UTC (10 years, 4 months ago) by jfenwick
File MIME type: text/plain
File size: 13209 byte(s)
Updating copyright notices
1
2 /*******************************************************
3 *
4 * Copyright (c) 2003-2009 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 if (scan_ret!=3)
95 {
96 MEMFREE(Ip);
97 MEMFREE(Jp);
98 MEMFREE(val);
99 fclose(f);
100 return -1;
101 }
102 /*FSCANF_CHECK(scan_ret, "fscanf: mm_read_unsymmetric_sparse");*/
103 Ip[i]--; /* adjust from 1-based to 0-based */
104 Jp[i]--;
105 }
106 fclose(f);
107
108 return 0;
109 }
110
111 int mm_is_valid(MM_typecode matcode)
112 {
113 if (!mm_is_matrix(matcode)) return 0;
114 if (mm_is_dense(matcode) && mm_is_pattern(matcode)) return 0;
115 if (mm_is_real(matcode) && mm_is_hermitian(matcode)) return 0;
116 if (mm_is_pattern(matcode) && (mm_is_hermitian(matcode) ||
117 mm_is_skew(matcode))) return 0;
118 return 1;
119 }
120
121 int mm_read_banner(FILE *f, MM_typecode *matcode)
122 {
123 char line[MM_MAX_LINE_LENGTH];
124 char banner[MM_MAX_TOKEN_LENGTH];
125 char mtx[MM_MAX_TOKEN_LENGTH];
126 char crd[MM_MAX_TOKEN_LENGTH];
127 char data_type[MM_MAX_TOKEN_LENGTH];
128 char storage_scheme[MM_MAX_TOKEN_LENGTH];
129 char *p;
130
131
132 mm_clear_typecode(matcode);
133
134 if (fgets(line, MM_MAX_LINE_LENGTH, f) == NULL)
135 return MM_PREMATURE_EOF;
136
137 if (sscanf(line, "%s %s %s %s %s", banner, mtx, crd, data_type,
138 storage_scheme) != 5)
139 return MM_PREMATURE_EOF;
140
141 /* convert to lower case */
142 for (p=mtx; *p!='\0'; *p= (char) tolower(*p),p++);
143 for (p=crd; *p!='\0'; *p= (char) tolower(*p),p++);
144 for (p=data_type; *p!='\0'; *p= (char) tolower(*p),p++);
145 for (p=storage_scheme; *p!='\0'; *p= (char) tolower(*p),p++);
146
147 /* check for banner */
148 if (strncmp(banner, MatrixMarketBanner, strlen(MatrixMarketBanner)) != 0)
149 return MM_NO_HEADER;
150
151 /* first field should be "mtx" */
152 if (strcmp(mtx, MM_MTX_STR) != 0)
153 return MM_UNSUPPORTED_TYPE;
154 mm_set_matrix(matcode);
155
156
157 /* second field describes whether this is a sparse matrix (in coordinate
158 storgae) or a dense array */
159
160
161 if (strcmp(crd, MM_SPARSE_STR) == 0)
162 mm_set_sparse(matcode);
163 else
164 if (strcmp(crd, MM_DENSE_STR) == 0)
165 mm_set_dense(matcode);
166 else
167 return MM_UNSUPPORTED_TYPE;
168
169
170 /* third field */
171
172 if (strcmp(data_type, MM_REAL_STR) == 0)
173 mm_set_real(matcode);
174 else
175 if (strcmp(data_type, MM_COMPLEX_STR) == 0)
176 mm_set_complex(matcode);
177 else
178 if (strcmp(data_type, MM_PATTERN_STR) == 0)
179 mm_set_pattern(matcode);
180 else
181 if (strcmp(data_type, MM_INT_STR) == 0)
182 mm_set_integer(matcode);
183 else
184 return MM_UNSUPPORTED_TYPE;
185
186
187 /* fourth field */
188
189 if (strcmp(storage_scheme, MM_GENERAL_STR) == 0)
190 mm_set_general(matcode);
191 else
192 if (strcmp(storage_scheme, MM_SYMM_STR) == 0)
193 mm_set_symmetric(matcode);
194 else
195 if (strcmp(storage_scheme, MM_HERM_STR) == 0)
196 mm_set_hermitian(matcode);
197 else
198 if (strcmp(storage_scheme, MM_SKEW_STR) == 0)
199 mm_set_skew(matcode);
200 else
201 return MM_UNSUPPORTED_TYPE;
202
203
204 return 0;
205 }
206
207 int mm_write_mtx_crd_size(FILE *f, int M, int N, int nz)
208 {
209 if (fprintf(f, "%d %d %d\n", M, N, nz) < 0)
210 return MM_COULD_NOT_WRITE_FILE;
211 else
212 return 0;
213 }
214
215 int mm_read_mtx_crd_size(FILE *f, int *M, int *N, int *nz )
216 {
217 char line[MM_MAX_LINE_LENGTH];
218 int num_items_read;
219
220 /* set return null parameter values, in case we exit with errors */
221 *M = *N = *nz = 0;
222
223 /* now continue scanning until you reach the end-of-comments */
224 do
225 {
226 if (fgets(line,MM_MAX_LINE_LENGTH,f) == NULL)
227 return MM_PREMATURE_EOF;
228 }while (line[0] == '%');
229
230 /* line[] is either blank or has M,N, nz */
231 if (sscanf(line, "%d %d %d", M, N, nz) == 3)
232 return 0;
233
234 else
235 do
236 {
237 num_items_read = fscanf(f, "%d %d %d", M, N, nz);
238 if (num_items_read == EOF) return MM_PREMATURE_EOF;
239 }
240 while (num_items_read != 3);
241
242 return 0;
243 }
244
245
246 int mm_read_mtx_array_size(FILE *f, int *M, int *N)
247 {
248 char line[MM_MAX_LINE_LENGTH];
249 int num_items_read;
250 /* set return null parameter values, in case we exit with errors */
251 *M = *N = 0;
252
253 /* now continue scanning until you reach the end-of-comments */
254 do
255 {
256 if (fgets(line,MM_MAX_LINE_LENGTH,f) == NULL)
257 return MM_PREMATURE_EOF;
258 }while (line[0] == '%');
259
260 /* line[] is either blank or has M,N, nz */
261 if (sscanf(line, "%d %d", M, N) == 2)
262 return 0;
263
264 else /* we have a blank line */
265 do
266 {
267 num_items_read = fscanf(f, "%d %d", M, N);
268 if (num_items_read == EOF) return MM_PREMATURE_EOF;
269 }
270 while (num_items_read != 2);
271
272 return 0;
273 }
274
275 int mm_write_mtx_array_size(FILE *f, int M, int N)
276 {
277 if (fprintf(f, "%d %d\n", M, N) < 0)
278 return MM_COULD_NOT_WRITE_FILE;
279 else
280 return 0;
281 }
282
283
284
285 /*-------------------------------------------------------------------------*/
286
287 /******************************************************************/
288 /* use when Ip[], Jp[], and val[] are already allocated */
289 /******************************************************************/
290
291 int mm_read_mtx_crd_data(FILE *f, int M, int N, int nz, int Ip[], int Jp[],
292 double val[], MM_typecode matcode)
293 {
294 int i;
295 if (mm_is_complex(matcode))
296 {
297 for (i=0; i<nz; i++)
298 if (fscanf(f, "%d %d %lg %lg", &Ip[i], &Jp[i], &val[2*i], &val[2*i+1])
299 != 4) return MM_PREMATURE_EOF;
300 }
301 else if (mm_is_real(matcode))
302 {
303 for (i=0; i<nz; i++)
304 {
305 if (fscanf(f, "%d %d %lg\n", &Ip[i], &Jp[i], &val[i])
306 != 3) return MM_PREMATURE_EOF;
307
308 }
309 }
310
311 else if (mm_is_pattern(matcode))
312 {
313 for (i=0; i<nz; i++)
314 if (fscanf(f, "%d %d", &Ip[i], &Jp[i])
315 != 2) return MM_PREMATURE_EOF;
316 }
317 else
318 return MM_UNSUPPORTED_TYPE;
319
320 return 0;
321
322 }
323
324 int mm_read_mtx_crd_entry(FILE *f, int *Ip, int *Jp,
325 double *real, double *imag, MM_typecode matcode)
326 {
327 if (mm_is_complex(matcode))
328 {
329 if (fscanf(f, "%d %d %lg %lg", Ip, Jp, real, imag)
330 != 4) return MM_PREMATURE_EOF;
331 }
332 else if (mm_is_real(matcode))
333 {
334 if (fscanf(f, "%d %d %lg\n", Ip, Jp, real)
335 != 3) return MM_PREMATURE_EOF;
336
337 }
338
339 else if (mm_is_pattern(matcode))
340 {
341 if (fscanf(f, "%d %d", Ip, Jp) != 2) return MM_PREMATURE_EOF;
342 }
343 else
344 return MM_UNSUPPORTED_TYPE;
345
346 return 0;
347
348 }
349
350
351 /************************************************************************
352 mm_read_mtx_crd() fills M, N, nz, array of values, and return
353 type code, e.g. 'MCRS'
354
355 if matrix is complex, values[] is of size 2*nz,
356 (nz pairs of real/imaginary values)
357 ************************************************************************/
358
359 int mm_read_mtx_crd(char *fname, int *M, int *N, int *nz, int **Ip, int **Jp,
360 double **val, MM_typecode *matcode)
361 {
362 int ret_code;
363 FILE *f;
364
365 if (strcmp(fname, "stdin") == 0) f=stdin;
366 else
367 if ((f = fopen(fname, "r")) == NULL)
368 return MM_COULD_NOT_READ_FILE;
369
370
371 if ((ret_code = mm_read_banner(f, matcode)) != 0)
372 return ret_code;
373
374 if (!(mm_is_valid(*matcode) && mm_is_sparse(*matcode) &&
375 mm_is_matrix(*matcode)))
376 return MM_UNSUPPORTED_TYPE;
377
378 if ((ret_code = mm_read_mtx_crd_size(f, M, N, nz)) != 0)
379 return ret_code;
380
381
382 *Ip = MEMALLOC(*nz, int);
383 *Jp = MEMALLOC(*nz, int);
384 *val = NULL;
385
386 if (mm_is_complex(*matcode))
387 {
388 *val = MEMALLOC(*nz * 2, double);
389 ret_code = mm_read_mtx_crd_data(f, *M, *N, *nz, *Ip, *Jp, *val,
390 *matcode);
391 if (ret_code != 0) return ret_code;
392 }
393 else if (mm_is_real(*matcode))
394 {
395 *val = MEMALLOC(*nz, double);
396 ret_code = mm_read_mtx_crd_data(f, *M, *N, *nz, *Ip, *Jp, *val,
397 *matcode);
398 if (ret_code != 0) return ret_code;
399 }
400
401 else if (mm_is_pattern(*matcode))
402 {
403 ret_code = mm_read_mtx_crd_data(f, *M, *N, *nz, *Ip, *Jp, *val,
404 *matcode);
405 if (ret_code != 0) return ret_code;
406 }
407
408 if (f != stdin) fclose(f);
409 return 0;
410 }
411
412 int mm_write_banner(FILE *f, MM_typecode matcode)
413 {
414 int ret_code = fprintf(f, "%s %s\n", MatrixMarketBanner, mm_typecode_to_str(matcode));
415 if (ret_code < 0 )
416 return MM_COULD_NOT_WRITE_FILE;
417 else
418 return 0;
419 }
420
421 int mm_write_mtx_crd(char fname[], int M, int N, int nz, int Ip[], int Jp[],
422 double val[], MM_typecode matcode)
423 {
424 FILE *f;
425 int i;
426
427 if (strcmp(fname, "stdout") == 0)
428 f = stdout;
429 else
430 if ((f = fopen(fname, "w")) == NULL)
431 return MM_COULD_NOT_WRITE_FILE;
432
433 /* print banner followed by typecode */
434 fprintf(f, "%s ", MatrixMarketBanner);
435 fprintf(f, "%s\n", mm_typecode_to_str(matcode));
436
437 /* print matrix sizes and nonzeros */
438 fprintf(f, "%d %d %d\n", M, N, nz);
439
440 /* print values */
441 if (mm_is_pattern(matcode))
442 for (i=0; i<nz; i++)
443 fprintf(f, "%d %d\n", Ip[i], Jp[i]);
444 else
445 if (mm_is_real(matcode))
446 for (i=0; i<nz; i++)
447 fprintf(f, "%d %d %20.16g\n", Ip[i], Jp[i], val[i]);
448 else
449 if (mm_is_complex(matcode))
450 for (i=0; i<nz; i++)
451 fprintf(f, "%d %d %20.16g %20.16g\n", Ip[i], Jp[i], val[2*i],
452 val[2*i+1]);
453 else
454 {
455 if (f != stdout) fclose(f);
456 return MM_UNSUPPORTED_TYPE;
457 }
458
459 if (f !=stdout) fclose(f);
460
461 return 0;
462 }
463
464
465 char *mm_typecode_to_str(MM_typecode matcode)
466 {
467 static char buffer[MM_MAX_LINE_LENGTH];
468 char *types[4];
469
470 /* check for MTX type */
471 if (mm_is_matrix(matcode))
472 types[0] = MM_MTX_STR;
473 else
474 return NULL;
475
476 /* check for CRD or ARR matrix */
477 if (mm_is_sparse(matcode))
478 types[1] = MM_SPARSE_STR;
479 else
480 if (mm_is_dense(matcode))
481 types[1] = MM_DENSE_STR;
482 else
483 return NULL;
484
485 /* check for element data type */
486 if (mm_is_real(matcode))
487 types[2] = MM_REAL_STR;
488 else
489 if (mm_is_complex(matcode))
490 types[2] = MM_COMPLEX_STR;
491 else
492 if (mm_is_pattern(matcode))
493 types[2] = MM_PATTERN_STR;
494 else
495 if (mm_is_integer(matcode))
496 types[2] = MM_INT_STR;
497 else
498 return NULL;
499
500
501 /* check for symmetry type */
502 if (mm_is_general(matcode))
503 types[3] = MM_GENERAL_STR;
504 else
505 if (mm_is_symmetric(matcode))
506 types[3] = MM_SYMM_STR;
507 else
508 if (mm_is_hermitian(matcode))
509 types[3] = MM_HERM_STR;
510 else
511 if (mm_is_skew(matcode))
512 types[3] = MM_SKEW_STR;
513 else
514 return NULL;
515
516 sprintf(buffer,"%s %s %s %s", types[0], types[1], types[2], types[3]);
517 return & buffer[0];
518
519 }
520
521 /*
522 * $Log$
523 * Revision 1.1 2004/10/26 06:53:59 jgs
524 * Initial revision
525 *
526 * Revision 1.1 2004/07/02 00:48:35 gross
527 * matrix market io function added
528 *
529 *
530 */

Properties

Name Value
svn:eol-style native
svn:keywords Author Date Id Revision

  ViewVC Help
Powered by ViewVC 1.1.26