Actual source code: mmio.c

  1: /*
  2:    Matrix Market I/O library for ANSI C

  4:    See https://math.nist.gov/MatrixMarket/ for details.
  5: */

  7: #include <stdlib.h>
  8: #include <stdio.h>
  9: #include <string.h>
 10: #include <ctype.h>

 12: #include "mmio.h"

 14: static char mm_buffer[MM_MAX_LINE_LENGTH];

 16: int mm_read_unsymmetric_sparse(const char *fname, int *M_, int *N_, int *nz_, double **val_, int **I_, int **J_)
 17: {
 18:   FILE       *f;
 19:   MM_typecode matcode;
 20:   int         M, N, nz;
 21:   int         i;
 22:   double     *val;
 23:   int        *ia, *ja;

 25:   if ((f = fopen(fname, "r")) == NULL) return -1;

 27:   if (mm_read_banner(f, &matcode) != 0) {
 28:     printf("mm_read_unsymmetric: Could not process Matrix Market banner ");
 29:     printf(" in file [%s]\n", fname);
 30:     return -1;
 31:   }

 33:   if (!(mm_is_real(matcode) && mm_is_matrix(matcode) && mm_is_sparse(matcode))) {
 34:     fprintf(stderr, "This application does not support ");
 35:     fprintf(stderr, "Market Market type: [%s]\n", mm_typecode_to_str(matcode));
 36:     return -1;
 37:   }

 39:   /* find out size of sparse matrix: M, N, nz .... */

 41:   if (mm_read_mtx_crd_size(f, &M, &N, &nz) != 0) {
 42:     fprintf(stderr, "read_unsymmetric_sparse(): could not parse matrix size.\n");
 43:     return -1;
 44:   }

 46:   *M_  = M;
 47:   *N_  = N;
 48:   *nz_ = nz;

 50:   /* reserve memory for matrices */

 52:   ia  = (int *)malloc(nz * sizeof(int));
 53:   ja  = (int *)malloc(nz * sizeof(int));
 54:   val = (double *)malloc(nz * sizeof(double));

 56:   *val_ = val;
 57:   *I_   = ia;
 58:   *J_   = ja;

 60:   /* NOTE: when reading in doubles, ANSI C requires the use of the "l"  */
 61:   /*   specifier as in "%lg", "%lf", "%le", otherwise errors will occur */
 62:   /*  (ANSI C X3.159-1989, Sec. 4.9.6.2, p. 136 lines 13-15)            */

 64:   for (i = 0; i < nz; i++) {
 65:     if (fscanf(f, "%d %d %lg\n", &ia[i], &ja[i], &val[i]) != 3) {
 66:       fprintf(stderr, "read_unsymmetric_sparse(): could not parse i, j and nonzero.\n");
 67:       return -1;
 68:     }
 69:     ia[i]--; /* adjust from 1-based to 0-based */
 70:     ja[i]--;
 71:   }
 72:   fclose(f);

 74:   return 0;
 75: }

 77: int mm_is_valid(MM_typecode matcode)
 78: {
 79:   if (!mm_is_matrix(matcode)) return 0;
 80:   if (mm_is_dense(matcode) && mm_is_pattern(matcode)) return 0;
 81:   if (mm_is_real(matcode) && mm_is_hermitian(matcode)) return 0;
 82:   if (mm_is_pattern(matcode) && (mm_is_hermitian(matcode) || mm_is_skew(matcode))) return 0;
 83:   return 1;
 84: }

 86: int mm_read_banner(FILE *f, MM_typecode *matcode)
 87: {
 88:   char  line[MM_MAX_LINE_LENGTH];
 89:   char  banner[MM_MAX_TOKEN_LENGTH];
 90:   char  mtx[MM_MAX_TOKEN_LENGTH];
 91:   char  crd[MM_MAX_TOKEN_LENGTH];
 92:   char  data_type[MM_MAX_TOKEN_LENGTH];
 93:   char  storage_scheme[MM_MAX_TOKEN_LENGTH];
 94:   char *p;

 96:   mm_clear_typecode(matcode);

 98:   if (fgets(line, MM_MAX_LINE_LENGTH, f) == NULL) return MM_PREMATURE_EOF;

100:   if (sscanf(line, "%s %s %s %s %s", banner, mtx, crd, data_type, storage_scheme) != 5) return MM_PREMATURE_EOF;

102:   for (p = mtx; *p != '\0'; *p = (char)tolower(*p), p++); /* convert to lower case */
103:   for (p = crd; *p != '\0'; *p = (char)tolower(*p), p++);
104:   for (p = data_type; *p != '\0'; *p = (char)tolower(*p), p++);
105:   for (p = storage_scheme; *p != '\0'; *p = (char)tolower(*p), p++);

107:   /* check for banner */
108:   if (strncmp(banner, MatrixMarketBanner, strlen(MatrixMarketBanner)) != 0) return MM_NO_HEADER;

110:   /* first field should be "mtx" */
111:   if (strcmp(mtx, MM_MTX_STR) != 0) return MM_UNSUPPORTED_TYPE;
112:   mm_set_matrix(matcode);

114:   /* second field describes whether this is a sparse matrix (in coordinate
115:             storgae) or a dense array */

117:   if (strcmp(crd, MM_SPARSE_STR) == 0) mm_set_sparse(matcode);
118:   else if (strcmp(crd, MM_DENSE_STR) == 0) mm_set_dense(matcode);
119:   else return MM_UNSUPPORTED_TYPE;

121:   /* third field */

123:   if (strcmp(data_type, MM_REAL_STR) == 0) mm_set_real(matcode);
124:   else if (strcmp(data_type, MM_COMPLEX_STR) == 0) mm_set_complex(matcode);
125:   else if (strcmp(data_type, MM_PATTERN_STR) == 0) mm_set_pattern(matcode);
126:   else if (strcmp(data_type, MM_INT_STR) == 0) mm_set_integer(matcode);
127:   else return MM_UNSUPPORTED_TYPE;

129:   /* fourth field */

131:   if (strcmp(storage_scheme, MM_GENERAL_STR) == 0) mm_set_general(matcode);
132:   else if (strcmp(storage_scheme, MM_SYMM_STR) == 0) mm_set_symmetric(matcode);
133:   else if (strcmp(storage_scheme, MM_HERM_STR) == 0) mm_set_hermitian(matcode);
134:   else if (strcmp(storage_scheme, MM_SKEW_STR) == 0) mm_set_skew(matcode);
135:   else return MM_UNSUPPORTED_TYPE;

137:   return 0;
138: }

140: int mm_write_mtx_crd_size(FILE *f, int M, int N, int nz)
141: {
142:   if (fprintf(f, "%d %d %d\n", M, N, nz) < 0) return MM_COULD_NOT_WRITE_FILE;
143:   else return 0;
144: }

146: int mm_read_mtx_crd_size(FILE *f, int *M, int *N, int *nz)
147: {
148:   char line[MM_MAX_LINE_LENGTH];
149:   int  num_items_read;

151:   /* set return null parameter values, in case we exit with errors */
152:   *M = *N = *nz = 0;

154:   /* now continue scanning until you reach the end-of-comments */
155:   do {
156:     if (fgets(line, MM_MAX_LINE_LENGTH, f) == NULL) return MM_PREMATURE_EOF;
157:   } while (line[0] == '%');

159:   /* line[] is either blank or has M,N, nz */
160:   if (sscanf(line, "%d %d %d", M, N, nz) == 3) return 0;

162:   else do {
163:       num_items_read = fscanf(f, "%d %d %d", M, N, nz);
164:       if (num_items_read == EOF) return MM_PREMATURE_EOF;
165:     } while (num_items_read != 3);

167:   return 0;
168: }

170: int mm_read_mtx_array_size(FILE *f, int *M, int *N)
171: {
172:   char line[MM_MAX_LINE_LENGTH];
173:   int  num_items_read;
174:   /* set return null parameter values, in case we exit with errors */
175:   *M = *N = 0;

177:   /* now continue scanning until you reach the end-of-comments */
178:   do {
179:     if (fgets(line, MM_MAX_LINE_LENGTH, f) == NULL) return MM_PREMATURE_EOF;
180:   } while (line[0] == '%');

182:   /* line[] is either blank or has M,N, nz */
183:   if (sscanf(line, "%d %d", M, N) == 2) return 0;

185:   else /* we have a blank line */ do {
186:       num_items_read = fscanf(f, "%d %d", M, N);
187:       if (num_items_read == EOF) return MM_PREMATURE_EOF;
188:     } while (num_items_read != 2);

190:   return 0;
191: }

193: int mm_write_mtx_array_size(FILE *f, int M, int N)
194: {
195:   if (fprintf(f, "%d %d\n", M, N) < 0) return MM_COULD_NOT_WRITE_FILE;
196:   else return 0;
197: }

199: /*-------------------------------------------------------------------------*/

201: /******************************************************************/
202: /* use when ia[], ja[], and val[] are already allocated */
203: /******************************************************************/

205: int mm_read_mtx_crd_data(FILE *f, int M, int N, int nz, int ia[], int ja[], double val[], MM_typecode matcode)
206: {
207:   int i;
208:   if (mm_is_complex(matcode)) {
209:     for (i = 0; i < nz; i++)
210:       if (fscanf(f, "%d %d %lg %lg", &ia[i], &ja[i], &val[2 * i], &val[2 * i + 1]) != 4) return MM_PREMATURE_EOF;
211:   } else if (mm_is_real(matcode)) {
212:     for (i = 0; i < nz; i++) {
213:       if (fscanf(f, "%d %d %lg\n", &ia[i], &ja[i], &val[i]) != 3) return MM_PREMATURE_EOF;
214:     }
215:   }

217:   else if (mm_is_pattern(matcode)) {
218:     for (i = 0; i < nz; i++)
219:       if (fscanf(f, "%d %d", &ia[i], &ja[i]) != 2) return MM_PREMATURE_EOF;
220:   } else return MM_UNSUPPORTED_TYPE;

222:   return 0;
223: }

225: int mm_read_mtx_crd_entry(FILE *f, int *ia, int *ja, double *real, double *imag, MM_typecode matcode)
226: {
227:   if (mm_is_complex(matcode)) {
228:     if (fscanf(f, "%d %d %lg %lg", ia, ja, real, imag) != 4) return MM_PREMATURE_EOF;
229:   } else if (mm_is_real(matcode)) {
230:     if (fscanf(f, "%d %d %lg\n", ia, ja, real) != 3) return MM_PREMATURE_EOF;

232:   }

234:   else if (mm_is_pattern(matcode)) {
235:     if (fscanf(f, "%d %d", ia, ja) != 2) return MM_PREMATURE_EOF;
236:   } else return MM_UNSUPPORTED_TYPE;

238:   return 0;
239: }

241: /************************************************************************
242:     mm_read_mtx_crd()  fills M, N, nz, array of values, and return
243:                         type code, e.g. 'MCRS'

245:                         if matrix is complex, values[] is of size 2*nz,
246:                             (nz pairs of real/imaginary values)
247: ************************************************************************/

249: int mm_read_mtx_crd(char *fname, int *M, int *N, int *nz, int **ia, int **ja, double **val, MM_typecode *matcode)
250: {
251:   int   ret_code;
252:   FILE *f;

254:   if (strcmp(fname, "stdin") == 0) f = stdin;
255:   else if ((f = fopen(fname, "r")) == NULL) return MM_COULD_NOT_READ_FILE;

257:   if ((ret_code = mm_read_banner(f, matcode)) != 0) return ret_code;

259:   if (!(mm_is_valid(*matcode) && mm_is_sparse(*matcode) && mm_is_matrix(*matcode))) return MM_UNSUPPORTED_TYPE;

261:   if ((ret_code = mm_read_mtx_crd_size(f, M, N, nz)) != 0) return ret_code;

263:   *ia  = (int *)malloc(*nz * sizeof(int));
264:   *ja  = (int *)malloc(*nz * sizeof(int));
265:   *val = NULL;

267:   if (mm_is_complex(*matcode)) {
268:     *val     = (double *)malloc(*nz * 2 * sizeof(double));
269:     ret_code = mm_read_mtx_crd_data(f, *M, *N, *nz, *ia, *ja, *val, *matcode);
270:     if (ret_code != 0) return ret_code;
271:   } else if (mm_is_real(*matcode)) {
272:     *val     = (double *)malloc(*nz * sizeof(double));
273:     ret_code = mm_read_mtx_crd_data(f, *M, *N, *nz, *ia, *ja, *val, *matcode);
274:     if (ret_code != 0) return ret_code;
275:   }

277:   else if (mm_is_pattern(*matcode)) {
278:     ret_code = mm_read_mtx_crd_data(f, *M, *N, *nz, *ia, *ja, *val, *matcode);
279:     if (ret_code != 0) return ret_code;
280:   }

282:   if (f != stdin) fclose(f);
283:   return 0;
284: }

286: int mm_write_banner(FILE *f, MM_typecode matcode)
287: {
288:   char *str = mm_typecode_to_str(matcode);
289:   int   ret_code;

291:   ret_code = fprintf(f, "%s %s\n", MatrixMarketBanner, str);
292:   if (ret_code < 0) return MM_COULD_NOT_WRITE_FILE;
293:   else return 0;
294: }

296: int mm_write_mtx_crd(char fname[], int M, int N, int nz, int ia[], int ja[], double val[], MM_typecode matcode)
297: {
298:   FILE *f;
299:   int   i;

301:   if (strcmp(fname, "stdout") == 0) f = stdout;
302:   else if ((f = fopen(fname, "w")) == NULL) return MM_COULD_NOT_WRITE_FILE;

304:   /* print banner followed by typecode */
305:   fprintf(f, "%s ", MatrixMarketBanner);
306:   fprintf(f, "%s\n", mm_typecode_to_str(matcode));

308:   /* print matrix sizes and nonzeros */
309:   fprintf(f, "%d %d %d\n", M, N, nz);

311:   /* print values */
312:   if (mm_is_pattern(matcode))
313:     for (i = 0; i < nz; i++) fprintf(f, "%d %d\n", ia[i], ja[i]);
314:   else if (mm_is_real(matcode))
315:     for (i = 0; i < nz; i++) fprintf(f, "%d %d %20.16g\n", ia[i], ja[i], val[i]);
316:   else if (mm_is_complex(matcode))
317:     for (i = 0; i < nz; i++) fprintf(f, "%d %d %20.16g %20.16g\n", ia[i], ja[i], val[2 * i], val[2 * i + 1]);
318:   else {
319:     if (f != stdout) fclose(f);
320:     return MM_UNSUPPORTED_TYPE;
321:   }

323:   if (f != stdout) fclose(f);

325:   return 0;
326: }

328: char *mm_typecode_to_str(MM_typecode matcode)
329: {
330:   const char *types[4];

332:   /* check for MTX type */
333:   if (mm_is_matrix(matcode)) types[0] = MM_MTX_STR;
334:   else return NULL;

336:   /* check for CRD or ARR matrix */
337:   if (mm_is_sparse(matcode)) types[1] = MM_SPARSE_STR;
338:   else if (mm_is_dense(matcode)) types[1] = MM_DENSE_STR;
339:   else return NULL;

341:   /* check for element data type */
342:   if (mm_is_real(matcode)) types[2] = MM_REAL_STR;
343:   else if (mm_is_complex(matcode)) types[2] = MM_COMPLEX_STR;
344:   else if (mm_is_pattern(matcode)) types[2] = MM_PATTERN_STR;
345:   else if (mm_is_integer(matcode)) types[2] = MM_INT_STR;
346:   else return NULL;

348:   /* check for symmetry type */
349:   if (mm_is_general(matcode)) types[3] = MM_GENERAL_STR;
350:   else if (mm_is_symmetric(matcode)) types[3] = MM_SYMM_STR;
351:   else if (mm_is_hermitian(matcode)) types[3] = MM_HERM_STR;
352:   else if (mm_is_skew(matcode)) types[3] = MM_SKEW_STR;
353:   else return NULL;

355:   mm_buffer[0] = 0;
356:   snprintf(mm_buffer, sizeof(mm_buffer) / sizeof(mm_buffer[0]), "%s %s %s %s", types[0], types[1], types[2], types[3]);
357:   return mm_buffer;
358: }