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,ierr,size,*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;

 40:   PetscInitialize(&argc,&argv,(char *)0,help);
 41:   MPI_Comm_size(PETSC_COMM_WORLD,&size);

 44:   PetscOptionsBegin(PETSC_COMM_WORLD,NULL,"Matrix Market example options","");
 45:   {
 46:     PetscOptionsString("-fin","Input Matrix Market file","",filein,filein,sizeof(filein),&flag);
 48:     PetscOptionsString("-fout","Output file in petsc sparse binary format","",fileout,fileout,sizeof(fileout),&flag);
 50:     PetscOptionsBool("-aij_only","Use MATAIJ for all cases","",aijonly,&aijonly,NULL);
 51:     PetscOptionsFList("-permute","Permute matrix and vector to solving in new ordering","",MatOrderingList,ordering,ordering,sizeof(ordering),&permute);
 52:   }
 53:   PetscOptionsEnd();

 55:   /* Read in matrix */
 56:   PetscFOpen(PETSC_COMM_SELF,filein,"r",&file);


 60:   /*  This is how one can screen matrix types if their application */
 61:   /*  only supports a subset of the Matrix Market data types.      */

 64:   if (mm_is_symmetric(matcode)) symmetric = PETSC_TRUE;
 65:   if (mm_is_skew(matcode)) skew = PETSC_TRUE;
 66:   if (mm_is_real(matcode)) real = PETSC_TRUE;
 67:   if (mm_is_pattern(matcode)) pattern = PETSC_TRUE;

 69:   /* Find out size of sparse matrix .... */

 72:   mm_write_banner(stdout, matcode);
 73:   PetscPrintf(PETSC_COMM_SELF,"M: %d, N: %d, nnz: %d\n",M,N,nz);

 75:   /* Reseve memory for matrices */
 76:   PetscMalloc4(nz,&ia,nz,&ja,nz,&val,M,&rownz);
 77:   for (i=0; i<M; i++) rownz[i] = 1; /* Since we will add 0.0 to diagonal entries */

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

 83:   for (i=0; i<nz; i++) {
 84:     if (pattern) {
 85:       ninput = fscanf(file, "%d %d\n", &ia[i], &ja[i]);
 87:       val[i] = 1.0;
 88:     } else if (real) {
 89:       ninput = fscanf(file, "%d %d %lg\n", &ia[i], &ja[i], &val[i]);
 91:     }
 92:     ia[i]--; ja[i]--;     /* adjust from 1-based to 0-based */
 93:     if (ia[i] != ja[i]) { /* already counted the diagonals above */
 94:       if ((symmetric && aijonly) || skew) { /* transpose */
 95:         rownz[ia[i]]++;
 96:         rownz[ja[i]]++;
 97:       } else rownz[ia[i]]++;
 98:     }
 99:   }
100:   PetscFClose(PETSC_COMM_SELF,file);
101:   PetscPrintf(PETSC_COMM_SELF,"Reading matrix completes.\n");

103:   /* Create, preallocate, and then assemble the matrix */
104:   MatCreate(PETSC_COMM_SELF,&A);
105:   MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,M,N);

107:   if (symmetric && !aijonly) {
108:     MatSetType(A,MATSEQSBAIJ);
109:     MatSetFromOptions(A);
110:     MatSetUp(A);
111:     MatSeqSBAIJSetPreallocation(A,1,0,rownz);
112:     PetscObjectTypeCompare((PetscObject)A,MATSEQSBAIJ,&sametype);
114:   } else {
115:     MatSetType(A,MATSEQAIJ);
116:     MatSetFromOptions(A);
117:     MatSetUp(A);
118:     MatSeqAIJSetPreallocation(A,0,rownz);
119:     PetscObjectTypeCompare((PetscObject)A,MATSEQAIJ,&sametype);
121:   }

123:   /* Add zero to diagonals, in case the matrix missing diagonals */
124:   for (j=0; j<M; j++)  MatSetValues(A,1,&j,1,&j,&zero,INSERT_VALUES);
125:   /* Add values to the matrix, these correspond to lower triangular part for symmetric or skew matrices */
126:   for (j=0; j<nz; j++) MatSetValues(A,1,&ia[j],1,&ja[j],&val[j],INSERT_VALUES);

128:   /* Add values to upper triangular part for some cases */
129:   if (symmetric && aijonly) {
130:     /* MatrixMarket matrix stores symm matrix in lower triangular part. Take its transpose */
131:     for (j=0; j<nz; j++) MatSetValues(A,1,&ja[j],1,&ia[j],&val[j],INSERT_VALUES);
132:   }
133:   if (skew) {
134:     for (j=0; j<nz; j++) {
135:       val[j] = -val[j];
136:       MatSetValues(A,1,&ja[j],1,&ia[j],&val[j],INSERT_VALUES);
137:     }
138:   }

140:   MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);
141:   MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);

143:   if (permute) {
144:     Mat Aperm;
145:     MatGetOrdering(A,ordering,&rowperm,&colperm);
146:     MatPermute(A,rowperm,colperm,&Aperm);
147:     MatDestroy(&A);
148:     A    = Aperm;               /* Replace original operator with permuted version */
149:   }

151:   /* Write out matrix */
152:   PetscPrintf(PETSC_COMM_SELF,"Writing matrix to binary file %s using PETSc %s format ...\n",fileout,(symmetric && !aijonly)?"SBAIJ":"AIJ");
153:   PetscViewerBinaryOpen(PETSC_COMM_SELF,fileout,FILE_MODE_WRITE,&view);
154:   MatView(A,view);
155:   PetscViewerDestroy(&view);
156:   PetscPrintf(PETSC_COMM_SELF,"Writing matrix completes.\n");

158:   PetscFree4(ia,ja,val,rownz);
159:   MatDestroy(&A);
160:   ISDestroy(&rowperm);
161:   ISDestroy(&colperm);
162:   PetscFinalize();
163:   return 0;
164: }

166: /*TEST

168:    build:
169:       requires:  !complex double !defined(PETSC_USE_64BIT_INDICES)
170:       depends: ex72mmio.c

172:    test:
173:       suffix: 1
174:       args: -fin ${wPETSC_DIR}/share/petsc/datafiles/matrices/amesos2_test_mat0.mtx -fout petscmat.aij
175:       output_file: output/ex72_1.out

177:    test:
178:       suffix: 2
179:       args: -fin ${wPETSC_DIR}/share/petsc/datafiles/matrices/LFAT5.mtx -fout petscmat.sbaij
180:       output_file: output/ex72_2.out

182:    test:
183:       suffix: 3
184:       args: -fin ${wPETSC_DIR}/share/petsc/datafiles/matrices/m_05_05_crk.mtx -fout petscmat2.aij
185:       output_file: output/ex72_3.out

187:    test:
188:       suffix: 4
189:       args: -fin ${wPETSC_DIR}/share/petsc/datafiles/matrices/amesos2_test_mat0.mtx -fout petscmat.aij -permute rcm
190:       output_file: output/ex72_4.out
191: TEST*/