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 -aij_only allows to use MATAIJ for all cases.\n\\n"; 7 8 /* 9 * NOTES: 10 * 11 * 1) Matrix Market files are always 1-based, i.e. the index of the first 12 * element of a matrix is (1,1), not (0,0) as in C. ADJUST THESE 13 * OFFSETS ACCORDINGLY offsets accordingly when reading and writing 14 * to files. 15 * 16 * 2) ANSI C requires one to use the "l" format modifier when reading 17 * double precision floating point numbers in scanf() and 18 * its variants. For example, use "%lf", "%lg", or "%le" 19 * when reading doubles, otherwise errors will occur. 20 */ 21 #include <petscmat.h> 22 #include "ex72mmio.h" 23 24 int main(int argc,char **argv) 25 { 26 PetscInt ret_code; 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 PetscInt i,j,nz,ierr,size,*rownz; 34 PetscScalar *val,zero = 0.0; 35 PetscViewer view; 36 PetscBool sametype,flag,symmetric = PETSC_FALSE,skew = PETSC_FALSE,real = PETSC_FALSE,pattern = PETSC_FALSE,aijonly = PETSC_FALSE; 37 38 PetscInitialize(&argc,&argv,(char *)0,help); 39 ierr = MPI_Comm_size(PETSC_COMM_WORLD,&size);CHKERRQ(ierr); 40 if (size != 1) SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_WRONG_MPI_SIZE,"This is a uniprocessor example only!"); 41 42 ierr = PetscOptionsGetString(NULL,NULL,"-fin",filein,sizeof(filein),&flag);CHKERRQ(ierr); 43 if (!flag) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_USER_INPUT,"Please use -fin <filename> to specify the input file name!"); 44 ierr = PetscOptionsGetString(NULL,NULL,"-fout",fileout,sizeof(fileout),&flag);CHKERRQ(ierr); 45 if (!flag) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_USER_INPUT,"Please use -fout <filename> to specify the output file name!"); 46 ierr = PetscOptionsGetBool(NULL,NULL,"-aij_only",&aijonly,NULL);CHKERRQ(ierr); 47 48 /* Read in matrix */ 49 ierr = PetscFOpen(PETSC_COMM_SELF,filein,"r",&file);CHKERRQ(ierr); 50 51 if (mm_read_banner(file, &matcode) != 0) 52 SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Could not process Matrix Market banner."); 53 54 /* This is how one can screen matrix types if their application */ 55 /* only supports a subset of the Matrix Market data types. */ 56 if (!mm_is_matrix(matcode) || !mm_is_sparse(matcode)){ 57 SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Input must be a sparse matrix. Market Market type: [%s]\n", mm_typecode_to_str(matcode)); 58 } 59 60 if (mm_is_symmetric(matcode)) symmetric = PETSC_TRUE; 61 if (mm_is_skew(matcode)) skew = PETSC_TRUE; 62 if (mm_is_real(matcode)) real = PETSC_TRUE; 63 if (mm_is_pattern(matcode)) pattern = PETSC_TRUE; 64 65 /* Find out size of sparse matrix .... */ 66 if ((ret_code = mm_read_mtx_crd_size(file, &M, &N, &nz)) !=0) 67 SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Size of sparse matrix is wrong."); 68 69 ierr = mm_write_banner(stdout, matcode);CHKERRQ(ierr); 70 ierr = PetscPrintf(PETSC_COMM_SELF,"M: %d, N: %d, nnz: %d\n",M,N,nz);CHKERRQ(ierr); 71 72 /* Reseve memory for matrices */ 73 ierr = PetscMalloc4(nz,&ia,nz,&ja,nz,&val,M,&rownz);CHKERRQ(ierr); 74 for (i=0; i<M; i++) rownz[i] = 1; /* Since we will add 0.0 to diagonal entries */ 75 76 /* NOTE: when reading in doubles, ANSI C requires the use of the "l" */ 77 /* specifier as in "%lg", "%lf", "%le", otherwise errors will occur */ 78 /* (ANSI C X3.159-1989, Sec. 4.9.6.2, p. 136 lines 13-15) */ 79 80 for (i=0; i<nz; i++){ 81 if (pattern) { 82 ninput = fscanf(file, "%d %d\n", &ia[i], &ja[i]); 83 if (ninput < 2) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_FILE_UNEXPECTED,"Badly formatted input file\n"); 84 val[i] = 1.0; 85 } else if (real) { 86 ninput = fscanf(file, "%d %d %lg\n", &ia[i], &ja[i], &val[i]); 87 if (ninput < 3) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_FILE_UNEXPECTED,"Badly formatted input file\n"); 88 } 89 ia[i]--; ja[i]--; /* adjust from 1-based to 0-based */ 90 if (ia[i] != ja[i]) { /* already counted the diagonals above */ 91 if ((symmetric && aijonly) || skew) { /* transpose */ 92 rownz[ia[i]]++; 93 rownz[ja[i]]++; 94 } else rownz[ia[i]]++; 95 } 96 } 97 ierr = PetscFClose(PETSC_COMM_SELF,file);CHKERRQ(ierr); 98 ierr = PetscPrintf(PETSC_COMM_SELF,"Reading matrix completes.\n");CHKERRQ(ierr); 99 100 /* Create, preallocate, and then assemble the matrix */ 101 ierr = MatCreate(PETSC_COMM_SELF,&A);CHKERRQ(ierr); 102 ierr = MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,M,N);CHKERRQ(ierr); 103 104 if (symmetric && !aijonly) { 105 ierr = MatSetType(A,MATSEQSBAIJ);CHKERRQ(ierr); 106 ierr = MatSetFromOptions(A);CHKERRQ(ierr); 107 ierr = MatSetUp(A);CHKERRQ(ierr); 108 ierr = MatSeqSBAIJSetPreallocation(A,1,0,rownz);CHKERRQ(ierr); 109 ierr = PetscObjectTypeCompare((PetscObject)A,MATSEQSBAIJ,&sametype);CHKERRQ(ierr); 110 if (!sametype) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Only AIJ and SBAIJ are supported. Your mattype is not supported"); 111 } else { 112 ierr = MatSetType(A,MATSEQAIJ);CHKERRQ(ierr); 113 ierr = MatSetFromOptions(A);CHKERRQ(ierr); 114 ierr = MatSetUp(A);CHKERRQ(ierr); 115 ierr = MatSeqAIJSetPreallocation(A,0,rownz);CHKERRQ(ierr); 116 ierr = PetscObjectTypeCompare((PetscObject)A,MATSEQAIJ,&sametype);CHKERRQ(ierr); 117 if (!sametype) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Only AIJ and SBAIJ are supported. Your mattype is not supported"); 118 } 119 120 /* Add zero to diagonals, in case the matrix missing diagonals */ 121 for (j=0; j<M; j++) { 122 ierr = MatSetValues(A,1,&j,1,&j,&zero,INSERT_VALUES);CHKERRQ(ierr); 123 } 124 /* Add values to the matrix, these correspond to lower triangular part for symmetric or skew matrices */ 125 for (j=0; j<nz; j++) { 126 ierr = MatSetValues(A,1,&ia[j],1,&ja[j],&val[j],INSERT_VALUES);CHKERRQ(ierr); 127 } 128 129 /* Add values to upper triangular part for some cases */ 130 if (symmetric && aijonly) { 131 /* MatrixMarket matrix stores symm matrix in lower triangular part. Take its transpose */ 132 for (j=0; j<nz; j++) { 133 ierr = MatSetValues(A,1,&ja[j],1,&ia[j],&val[j],INSERT_VALUES);CHKERRQ(ierr); 134 } 135 } 136 if (skew) { 137 for (j=0; j<nz; j++) { 138 val[j] = -val[j]; 139 ierr = MatSetValues(A,1,&ja[j],1,&ia[j],&val[j],INSERT_VALUES);CHKERRQ(ierr); 140 } 141 } 142 143 ierr = MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 144 ierr = MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 145 146 /* Write out matrix */ 147 ierr = PetscPrintf(PETSC_COMM_SELF,"Writing matrix to binary file %s using PETSc %s format ...\n",fileout,(symmetric && !aijonly)?"SBAIJ":"AIJ");CHKERRQ(ierr); 148 ierr = PetscViewerBinaryOpen(PETSC_COMM_SELF,fileout,FILE_MODE_WRITE,&view);CHKERRQ(ierr); 149 ierr = MatView(A,view);CHKERRQ(ierr); 150 ierr = PetscViewerDestroy(&view);CHKERRQ(ierr); 151 ierr = PetscPrintf(PETSC_COMM_SELF,"Writing matrix completes.\n");CHKERRQ(ierr); 152 153 ierr = PetscFree4(ia,ja,val,rownz);CHKERRQ(ierr); 154 ierr = MatDestroy(&A);CHKERRQ(ierr); 155 ierr = PetscFinalize();CHKERRQ(ierr); 156 return 0; 157 } 158 159 /*TEST 160 161 build: 162 requires: !complex double !define(PETSC_USE_64BIT_INDICES) 163 depends: ex72mmio.c 164 165 test: 166 suffix: 1 167 args: -fin ${wPETSC_DIR}/share/petsc/datafiles/matrices/amesos2_test_mat0.mtx -fout petscmat.aij 168 output_file: output/ex72_1.out 169 170 test: 171 suffix: 2 172 args: -fin ${wPETSC_DIR}/share/petsc/datafiles/matrices/LFAT5.mtx -fout petscmat.sbaij 173 output_file: output/ex72_2.out 174 175 test: 176 suffix: 3 177 args: -fin ${wPETSC_DIR}/share/petsc/datafiles/matrices/m_05_05_crk.mtx -fout petscmat2.aij 178 output_file: output/ex72_3.out 179 TEST*/ 180