Graphviz 12.0.1~dev.20240716.0800
Loading...
Searching...
No Matches
matrix_market.c
Go to the documentation of this file.
1/*************************************************************************
2 * Copyright (c) 2011 AT&T Intellectual Property
3 * All rights reserved. This program and the accompanying materials
4 * are made available under the terms of the Eclipse Public License v1.0
5 * which accompanies this distribution, and is available at
6 * https://www.eclipse.org/legal/epl-v10.html
7 *
8 * Contributors: Details at https://graphviz.org
9 *************************************************************************/
10
11#include <sparse/SparseMatrix.h>
12#include "mmio.h"
13#include "matrix_market.h"
14#include <cgraph/alloc.h>
15#include <cgraph/unreachable.h>
16#include <assert.h>
17
18static int mm_get_type(MM_typecode typecode)
19{
20 if (mm_is_complex(typecode)) {
21 return MATRIX_TYPE_COMPLEX;
22 } else if (mm_is_real(typecode)) {
23 return MATRIX_TYPE_REAL;
24 } else if (mm_is_integer(typecode)) {
25 return MATRIX_TYPE_INTEGER;
26 } else if (mm_is_pattern(typecode)) {
27 return MATRIX_TYPE_PATTERN;
28 }
29 return MATRIX_TYPE_UNKNOWN;
30}
31
32SparseMatrix SparseMatrix_import_matrix_market(FILE * f)
33{
34 int ret_code, type;
35 MM_typecode matcode;
36 double *val = NULL, *v;
37 int *vali = NULL, i, m, n, *I = NULL, *J = NULL, nz;
38 void *vp = NULL;
39 SparseMatrix A = NULL;
40 int nzold;
41 int c;
42
43 if ((c = fgetc(f)) != '%') {
44 ungetc(c, f);
45 return NULL;
46 }
47 ungetc(c, f);
48 if (mm_read_banner(f, &matcode) != 0) {
49#ifdef DEBUG
50 printf("Could not process Matrix Market banner.\n");
51#endif
52 return NULL;
53 }
54
55
56 /* This is how one can screen matrix types if their application */
57 /* only supports a subset of the Matrix Market data types. */
58
59 if (!mm_is_matrix(matcode) || !mm_is_sparse(matcode)) {
60 UNREACHABLE();
61 }
62
63 /* find out size of sparse matrix .... */
64 if ((ret_code = mm_read_mtx_crd_size(f, &m, &n, &nz)) != 0) {
65 assert(0);
66 return NULL;
67 }
68 /* reseve memory for matrices */
69
70 I = gv_calloc(nz, sizeof(int));
71 J = gv_calloc(nz, sizeof(int));
72
73 /* NOTE: when reading in doubles, ANSI C requires the use of the "l" */
74 /* specifier as in "%lg", "%lf", "%le", otherwise errors will occur */
75 /* (ANSI C X3.159-1989, Sec. 4.9.6.2, p. 136 lines 13-15) */
76 type = mm_get_type(matcode);
77 switch (type) {
78 case MATRIX_TYPE_REAL:
79 val = gv_calloc(nz, sizeof(double));
80 for (i = 0; i < nz; i++) {
81 int num = fscanf(f, "%d %d %lg\n", &I[i], &J[i], &val[i]);
82 (void)num;
83 assert(num == 3);
84 I[i]--; /* adjust from 1-based to 0-based */
85 J[i]--;
86 }
87 if (mm_is_symmetric(matcode)) {
88 I = gv_recalloc(I, nz, 2 * nz, sizeof(int));
89 J = gv_recalloc(J, nz, 2 * nz, sizeof(int));
90 val = gv_recalloc(val, nz, 2 * nz, sizeof(double));
91 nzold = nz;
92 for (i = 0; i < nzold; i++) {
93 if (I[i] != J[i]) {
94 I[nz] = J[i];
95 J[nz] = I[i];
96 val[nz++] = val[i];
97 }
98 }
99 } else if (mm_is_skew(matcode)) {
100 I = gv_recalloc(I, nz, 2 * nz, sizeof(int));
101 J = gv_recalloc(J, nz, 2 * nz, sizeof(int));
102 val = gv_recalloc(val, nz, 2 * nz, sizeof(double));
103 nzold = nz;
104 for (i = 0; i < nzold; i++) {
105 assert(I[i] != J[i]); /* skew symm has no diag */
106 I[nz] = J[i];
107 J[nz] = I[i];
108 val[nz++] = -val[i];
109 }
110 } else {
111 assert(!mm_is_hermitian(matcode));
112 }
113 vp = val;
114 break;
115 case MATRIX_TYPE_INTEGER:
116 vali = gv_calloc(nz, sizeof(int));
117 for (i = 0; i < nz; i++) {
118 int num = fscanf(f, "%d %d %d\n", &I[i], &J[i], &vali[i]);
119 (void)num;
120 assert(num == 3);
121 I[i]--; /* adjust from 1-based to 0-based */
122 J[i]--;
123 }
124 if (mm_is_symmetric(matcode)) {
125 I = gv_recalloc(I, nz, 2 * nz, sizeof(int));
126 J = gv_recalloc(J, nz, 2 * nz, sizeof(int));
127 vali = gv_recalloc(vali, nz, 2 * nz, sizeof(int));
128 nzold = nz;
129 for (i = 0; i < nzold; i++) {
130 if (I[i] != J[i]) {
131 I[nz] = J[i];
132 J[nz] = I[i];
133 vali[nz++] = vali[i];
134 }
135 }
136 } else if (mm_is_skew(matcode)) {
137 I = gv_recalloc(I, nz, 2 * nz, sizeof(int));
138 J = gv_recalloc(J, nz, 2 * nz, sizeof(int));
139 vali = gv_recalloc(vali, nz, 2 * nz, sizeof(int));
140 nzold = nz;
141 for (i = 0; i < nzold; i++) {
142 assert(I[i] != J[i]); /* skew symm has no diag */
143 I[nz] = J[i];
144 J[nz] = I[i];
145 vali[nz++] = -vali[i];
146 }
147 } else {
148 assert(!mm_is_hermitian(matcode));
149 }
150 vp = vali;
151 break;
152 case MATRIX_TYPE_PATTERN:
153 for (i = 0; i < nz; i++) {
154 int num = fscanf(f, "%d %d\n", &I[i], &J[i]);
155 (void)num;
156 assert(num == 2);
157 I[i]--; /* adjust from 1-based to 0-based */
158 J[i]--;
159 }
160 if (mm_is_symmetric(matcode) || mm_is_skew(matcode)) {
161 I = gv_recalloc(I, nz, 2 * nz, sizeof(int));
162 J = gv_recalloc(J, nz, 2 * nz, sizeof(int));
163 nzold = nz;
164 for (i = 0; i < nzold; i++) {
165 if (I[i] != J[i]) {
166 I[nz] = J[i];
167 J[nz++] = I[i];
168 }
169 }
170 } else {
171 assert(!mm_is_hermitian(matcode));
172 }
173 break;
174 case MATRIX_TYPE_COMPLEX:
175 val = gv_calloc(2 * nz, sizeof(double));
176 v = val;
177 for (i = 0; i < nz; i++) {
178 int num = fscanf(f, "%d %d %lg %lg\n", &I[i], &J[i], &v[0], &v[1]);
179 (void)num;
180 assert(num == 4);
181 v += 2;
182 I[i]--; /* adjust from 1-based to 0-based */
183 J[i]--;
184 }
185 if (mm_is_symmetric(matcode)) {
186 I = gv_recalloc(I, nz, 2 * nz, sizeof(int));
187 J = gv_recalloc(J, nz, 2 * nz, sizeof(int));
188 val = gv_recalloc(val, 2 * nz, 4 * nz, sizeof(double));
189 nzold = nz;
190 for (i = 0; i < nzold; i++) {
191 if (I[i] != J[i]) {
192 I[nz] = J[i];
193 J[nz] = I[i];
194 val[2 * nz] = val[2 * i];
195 val[2 * nz + 1] = val[2 * i + 1];
196 nz++;
197 }
198 }
199 } else if (mm_is_skew(matcode)) {
200 I = gv_recalloc(I, nz, 2 * nz, sizeof(int));
201 J = gv_recalloc(J, nz, 2 * nz, sizeof(int));
202 val = gv_recalloc(val, 2 * nz, 4 * nz, sizeof(double));
203 nzold = nz;
204 for (i = 0; i < nzold; i++) {
205 assert(I[i] != J[i]); /* skew symm has no diag */
206 I[nz] = J[i];
207 J[nz] = I[i];
208 val[2 * nz] = -val[2 * i];
209 val[2 * nz + 1] = -val[2 * i + 1];
210 nz++;
211
212 }
213 } else if (mm_is_hermitian(matcode)) {
214 I = gv_recalloc(I, nz, 2 * nz, sizeof(int));
215 J = gv_recalloc(J, nz, 2 * nz, sizeof(int));
216 val = gv_recalloc(val, 2 * nz, 4 * nz, sizeof(double));
217 nzold = nz;
218 for (i = 0; i < nzold; i++) {
219 if (I[i] != J[i]) {
220 I[nz] = J[i];
221 J[nz] = I[i];
222 val[2 * nz] = val[2 * i];
223 val[2 * nz + 1] = -val[2 * i + 1];
224 nz++;
225 }
226 }
227 }
228 vp = val;
229 break;
230 default:
231 return 0;
232 }
233
234 A = SparseMatrix_from_coordinate_arrays(nz, m, n, I, J, vp,
235 type, sizeof(double));
236 free(I);
237 free(J);
238 free(val);
239
240 if (mm_is_symmetric(matcode)) {
241 SparseMatrix_set_symmetric(A);
242 SparseMatrix_set_pattern_symmetric(A);
243 }
244
245
246 return A;
247}
SparseMatrix_from_coordinate_arrays
SparseMatrix SparseMatrix_from_coordinate_arrays(int nz, int m, int n, int *irn, int *jcn, void *val0, int type, size_t sz)
Definition SparseMatrix.c:694
SparseMatrix_set_pattern_symmetric
#define SparseMatrix_set_pattern_symmetric(A)
Definition SparseMatrix.h:111
MATRIX_TYPE_REAL
@ MATRIX_TYPE_REAL
Definition SparseMatrix.h:45
MATRIX_TYPE_UNKNOWN
@ MATRIX_TYPE_UNKNOWN
Definition SparseMatrix.h:45
MATRIX_TYPE_PATTERN
@ MATRIX_TYPE_PATTERN
Definition SparseMatrix.h:45
MATRIX_TYPE_COMPLEX
@ MATRIX_TYPE_COMPLEX
Definition SparseMatrix.h:45
MATRIX_TYPE_INTEGER
@ MATRIX_TYPE_INTEGER
Definition SparseMatrix.h:45
SparseMatrix_set_symmetric
#define SparseMatrix_set_symmetric(A)
Definition SparseMatrix.h:110
alloc.h
Memory allocation wrappers that exit on failure.
gv_recalloc
static void * gv_recalloc(void *ptr, size_t old_nmemb, size_t new_nmemb, size_t size)
Definition alloc.h:73
gv_calloc
static void * gv_calloc(size_t nmemb, size_t size)
Definition alloc.h:26
type
expr procedure type
Definition exparse.y:211
A
#define A(n, t)
Definition expr.h:76
I
#define I
Definition expr.h:71
free
void free(void *)
NULL
node NULL
Definition grammar.y:149
SparseMatrix_import_matrix_market
SparseMatrix SparseMatrix_import_matrix_market(FILE *f)
Definition matrix_market.c:32
mm_get_type
static int mm_get_type(MM_typecode typecode)
Definition matrix_market.c:18
matrix_market.h
import Matrix Market
mm_read_mtx_crd_size
int mm_read_mtx_crd_size(FILE *f, int *M, int *N, int *nz)
Definition mmio.c:100
mm_read_banner
int mm_read_banner(FILE *f, MM_typecode *matcode)
Definition mmio.c:28
mmio.h
Matrix Market I/O API
mm_is_complex
#define mm_is_complex(typecode)
Definition mmio.h:43
MM_typecode
char MM_typecode[4]
Definition mmio.h:29
mm_is_matrix
#define mm_is_matrix(typecode)
Definition mmio.h:36
mm_is_hermitian
#define mm_is_hermitian(typecode)
Definition mmio.h:51
mm_is_real
#define mm_is_real(typecode)
Definition mmio.h:44
mm_is_skew
#define mm_is_skew(typecode)
Definition mmio.h:50
mm_is_pattern
#define mm_is_pattern(typecode)
Definition mmio.h:45
mm_is_integer
#define mm_is_integer(typecode)
Definition mmio.h:46
mm_is_symmetric
#define mm_is_symmetric(typecode)
Definition mmio.h:48
mm_is_sparse
#define mm_is_sparse(typecode)
Definition mmio.h:38
SparseMatrix_struct
Definition SparseMatrix.h:29
UNREACHABLE
#define UNREACHABLE()
Definition unreachable.h:30