Actual source code: ex72.c
1: static char help[] = "Read a non-complex sparse matrix from a Matrix Market (v. 2.0) file\n\
2: and write it to a file in petsc sparse binary format. If the matrix is symmetric, the binary file is in \n\
3: PETSc MATSBAIJ format, otherwise it is in MATAIJ format \n\
4: Usage: ./ex72 -fin <infile> -fout <outfile> \n\
5: (See https://math.nist.gov/MatrixMarket/ for details.)\n\
6: The option -permute <natural,rcm,nd,...> permutes the matrix using the ordering type.\n\
7: The option -aij_only allows to use MATAIJ for all cases.\n\\n";
9: /*
10: NOTES:
12: 1) Matrix Market files are always 1-based, i.e. the index of the first
13: element of a matrix is (1,1), not (0,0) as in C. ADJUST THESE
14: OFFSETS ACCORDINGLY offsets accordingly when reading and writing
15: to files.
17: 2) ANSI C requires one to use the "l" format modifier when reading
18: double precision floating point numbers in scanf() and
19: its variants. For example, use "%lf", "%lg", or "%le"
20: when reading doubles, otherwise errors will occur.
21: */
22: #include <petscmat.h>
23: #include "ex72mmio.h"
25: int main(int argc, char **argv)
26: {
27: MM_typecode matcode;
28: FILE *file;
29: PetscInt M, N, ninput;
30: PetscInt *ia, *ja;
31: Mat A;
32: char filein[PETSC_MAX_PATH_LEN], fileout[PETSC_MAX_PATH_LEN];
33: char ordering[256] = MATORDERINGRCM;
34: PetscInt i, j, nz, *rownz;
35: PetscScalar *val, zero = 0.0;
36: PetscViewer view;
37: PetscBool sametype, flag, symmetric = PETSC_FALSE, skew = PETSC_FALSE, real = PETSC_FALSE, pattern = PETSC_FALSE, aijonly = PETSC_FALSE, permute = PETSC_FALSE;
38: IS rowperm = NULL, colperm = NULL;
39: PetscMPIInt size;
41: PetscInitialize(&argc, &argv, (char *)0, help);
42: MPI_Comm_size(PETSC_COMM_WORLD, &size);
45: PetscOptionsBegin(PETSC_COMM_WORLD, NULL, "Matrix Market example options", "");
46: {
47: PetscOptionsString("-fin", "Input Matrix Market file", "", filein, filein, sizeof(filein), &flag);
49: PetscOptionsString("-fout", "Output file in petsc sparse binary format", "", fileout, fileout, sizeof(fileout), &flag);
51: PetscOptionsBool("-aij_only", "Use MATAIJ for all cases", "", aijonly, &aijonly, NULL);
52: PetscOptionsFList("-permute", "Permute matrix and vector to solving in new ordering", "", MatOrderingList, ordering, ordering, sizeof(ordering), &permute);
53: }
54: PetscOptionsEnd();
56: /* Read in matrix */
57: PetscFOpen(PETSC_COMM_SELF, filein, "r", &file);
61: /* This is how one can screen matrix types if their application */
62: /* only supports a subset of the Matrix Market data types. */
65: if (mm_is_symmetric(matcode)) symmetric = PETSC_TRUE;
66: if (mm_is_skew(matcode)) skew = PETSC_TRUE;
67: if (mm_is_real(matcode)) real = PETSC_TRUE;
68: if (mm_is_pattern(matcode)) pattern = PETSC_TRUE;
70: /* Find out size of sparse matrix .... */
73: mm_write_banner(stdout, matcode);
74: PetscPrintf(PETSC_COMM_SELF, "M: %d, N: %d, nnz: %d\n", M, N, nz);
76: /* Reseve memory for matrices */
77: PetscMalloc4(nz, &ia, nz, &ja, nz, &val, M, &rownz);
78: for (i = 0; i < M; i++) rownz[i] = 1; /* Since we will add 0.0 to diagonal entries */
80: /* NOTE: when reading in doubles, ANSI C requires the use of the "l" */
81: /* specifier as in "%lg", "%lf", "%le", otherwise errors will occur */
82: /* (ANSI C X3.159-1989, Sec. 4.9.6.2, p. 136 lines 13-15) */
84: for (i = 0; i < nz; i++) {
85: if (pattern) {
86: ninput = fscanf(file, "%d %d\n", &ia[i], &ja[i]);
88: val[i] = 1.0;
89: } else if (real) {
90: ninput = fscanf(file, "%d %d %lg\n", &ia[i], &ja[i], &val[i]);
92: }
93: ia[i]--;
94: ja[i]--; /* adjust from 1-based to 0-based */
95: if (ia[i] != ja[i]) { /* already counted the diagonals above */
96: if ((symmetric && aijonly) || skew) { /* transpose */
97: rownz[ia[i]]++;
98: rownz[ja[i]]++;
99: } else rownz[ia[i]]++;
100: }
101: }
102: PetscFClose(PETSC_COMM_SELF, file);
103: PetscPrintf(PETSC_COMM_SELF, "Reading matrix completes.\n");
105: /* Create, preallocate, and then assemble the matrix */
106: MatCreate(PETSC_COMM_SELF, &A);
107: MatSetSizes(A, PETSC_DECIDE, PETSC_DECIDE, M, N);
109: if (symmetric && !aijonly) {
110: MatSetType(A, MATSEQSBAIJ);
111: MatSetFromOptions(A);
112: MatSetUp(A);
113: MatSeqSBAIJSetPreallocation(A, 1, 0, rownz);
114: PetscObjectTypeCompare((PetscObject)A, MATSEQSBAIJ, &sametype);
116: } else {
117: MatSetType(A, MATSEQAIJ);
118: MatSetFromOptions(A);
119: MatSetUp(A);
120: MatSeqAIJSetPreallocation(A, 0, rownz);
121: PetscObjectTypeCompare((PetscObject)A, MATSEQAIJ, &sametype);
123: }
125: /* Add zero to diagonals, in case the matrix missing diagonals */
126: for (j = 0; j < M; j++) MatSetValues(A, 1, &j, 1, &j, &zero, INSERT_VALUES);
127: /* Add values to the matrix, these correspond to lower triangular part for symmetric or skew matrices */
128: for (j = 0; j < nz; j++) MatSetValues(A, 1, &ia[j], 1, &ja[j], &val[j], INSERT_VALUES);
130: /* Add values to upper triangular part for some cases */
131: if (symmetric && aijonly) {
132: /* MatrixMarket matrix stores symm matrix in lower triangular part. Take its transpose */
133: for (j = 0; j < nz; j++) MatSetValues(A, 1, &ja[j], 1, &ia[j], &val[j], INSERT_VALUES);
134: }
135: if (skew) {
136: for (j = 0; j < nz; j++) {
137: val[j] = -val[j];
138: MatSetValues(A, 1, &ja[j], 1, &ia[j], &val[j], INSERT_VALUES);
139: }
140: }
142: MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY);
143: MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY);
145: if (permute) {
146: Mat Aperm;
147: MatGetOrdering(A, ordering, &rowperm, &colperm);
148: MatPermute(A, rowperm, colperm, &Aperm);
149: MatDestroy(&A);
150: A = Aperm; /* Replace original operator with permuted version */
151: }
153: /* Write out matrix */
154: PetscPrintf(PETSC_COMM_SELF, "Writing matrix to binary file %s using PETSc %s format ...\n", fileout, (symmetric && !aijonly) ? "SBAIJ" : "AIJ");
155: PetscViewerBinaryOpen(PETSC_COMM_SELF, fileout, FILE_MODE_WRITE, &view);
156: MatView(A, view);
157: PetscViewerDestroy(&view);
158: PetscPrintf(PETSC_COMM_SELF, "Writing matrix completes.\n");
160: PetscFree4(ia, ja, val, rownz);
161: MatDestroy(&A);
162: ISDestroy(&rowperm);
163: ISDestroy(&colperm);
164: PetscFinalize();
165: return 0;
166: }
168: /*TEST
170: build:
171: requires: !complex double !defined(PETSC_USE_64BIT_INDICES)
172: depends: ex72mmio.c
174: test:
175: suffix: 1
176: args: -fin ${wPETSC_DIR}/share/petsc/datafiles/matrices/amesos2_test_mat0.mtx -fout petscmat.aij
177: output_file: output/ex72_1.out
179: test:
180: suffix: 2
181: args: -fin ${wPETSC_DIR}/share/petsc/datafiles/matrices/LFAT5.mtx -fout petscmat.sbaij
182: output_file: output/ex72_2.out
184: test:
185: suffix: 3
186: args: -fin ${wPETSC_DIR}/share/petsc/datafiles/matrices/m_05_05_crk.mtx -fout petscmat2.aij
187: output_file: output/ex72_3.out
189: test:
190: suffix: 4
191: args: -fin ${wPETSC_DIR}/share/petsc/datafiles/matrices/amesos2_test_mat0.mtx -fout petscmat.aij -permute rcm
192: output_file: output/ex72_4.out
193: TEST*/