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