1 2 static char help[] = "Reads in a matrix in ASCII MATLAB format (I,J,A), read in vectors rhs and exact_solu in ASCII format.\n\ 3 Writes them using the PETSc sparse format.\n\ 4 Note: I and J start at 1, not 0, use -noshift if indices in file start with zero!\n\ 5 Input parameters are:\n\ 6 -Ain <filename> : input matrix in ascii format\n\ 7 -rhs <filename> : input rhs in ascii format\n\ 8 -solu <filename> : input true solution in ascii format\n\\n"; 9 10 /* 11 Example: ./ex78 -Ain Ain -rhs rhs -solu solu -noshift -mat_view 12 with the datafiles in the following format: 13 Ain (I and J start at 0): 14 ------------------------ 15 3 3 6 16 0 0 1.0 17 0 1 2.0 18 1 0 3.0 19 1 1 4.0 20 1 2 5.0 21 2 2 6.0 22 23 rhs 24 --- 25 0 3.0 26 1 12.0 27 2 6.0 28 29 solu 30 ---- 31 0 1.0 32 0 1.0 33 0 1.0 34 */ 35 36 #include <petscmat.h> 37 38 int main(int argc,char **args) 39 { 40 Mat A = NULL; 41 Vec b,u = NULL,u_tmp; 42 char Ain[PETSC_MAX_PATH_LEN],rhs[PETSC_MAX_PATH_LEN],solu[PETSC_MAX_PATH_LEN]; 43 int m,n = 0,nz,dummy; /* these are fscaned so kept as int */ 44 PetscInt i,col,row,shift = 1,sizes[3],nsizes; 45 PetscScalar val; 46 PetscReal res_norm; 47 FILE *Afile,*bfile,*ufile; 48 PetscViewer view; 49 PetscBool flg_A,flg_b,flg_u,flg; 50 PetscMPIInt size; 51 52 CHKERRQ(PetscInitialize(&argc,&args,(char*)0,help)); 53 CHKERRMPI(MPI_Comm_size(PETSC_COMM_WORLD,&size)); 54 PetscCheckFalse(size != 1,PETSC_COMM_WORLD,PETSC_ERR_SUP,"This is a uniprocessor example only!"); 55 56 /* Read in matrix, rhs and exact solution from ascii files */ 57 CHKERRQ(PetscOptionsGetString(NULL,NULL,"-Ain",Ain,sizeof(Ain),&flg_A)); 58 CHKERRQ(PetscOptionsHasName(NULL,NULL,"-noshift",&flg)); 59 if (flg) shift = 0; 60 if (flg_A) { 61 CHKERRQ(PetscPrintf(PETSC_COMM_SELF,"\n Read matrix in ascii format ...\n")); 62 CHKERRQ(PetscFOpen(PETSC_COMM_SELF,Ain,"r",&Afile)); 63 nsizes = 3; 64 CHKERRQ(PetscOptionsGetIntArray(NULL,NULL,"-nosizesinfile",sizes,&nsizes,&flg)); 65 if (flg) { 66 PetscCheckFalse(nsizes != 3,PETSC_COMM_WORLD,PETSC_ERR_USER,"Must pass in three m,n,nz as arguments for -nosizesinfile"); 67 m = sizes[0]; 68 n = sizes[1]; 69 nz = sizes[2]; 70 } else { 71 PetscCheckFalse(fscanf(Afile,"%d %d %d\n",&m,&n,&nz) != 3,PETSC_COMM_SELF,PETSC_ERR_FILE_UNEXPECTED,"Badly formatted input file"); 72 } 73 CHKERRQ(PetscPrintf(PETSC_COMM_SELF,"m: %d, n: %d, nz: %d \n", m,n,nz)); 74 PetscCheckFalse(m != n,PETSC_COMM_SELF,PETSC_ERR_ARG_SIZ, "Number of rows, cols must be same for this example"); 75 CHKERRQ(MatCreate(PETSC_COMM_SELF,&A)); 76 CHKERRQ(MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,m,n)); 77 CHKERRQ(MatSetFromOptions(A)); 78 CHKERRQ(MatSeqAIJSetPreallocation(A,nz/m,NULL)); 79 CHKERRQ(MatSetOption(A,MAT_NEW_NONZERO_ALLOCATION_ERR,PETSC_FALSE)); 80 81 for (i=0; i<nz; i++) { 82 PetscCheckFalse(fscanf(Afile,"%d %d %le\n",&row,&col,(double*)&val) != 3,PETSC_COMM_SELF,PETSC_ERR_FILE_UNEXPECTED,"Badly formatted input file"); 83 row -= shift; col -= shift; /* set index set starts at 0 */ 84 CHKERRQ(MatSetValues(A,1,&row,1,&col,&val,INSERT_VALUES)); 85 } 86 CHKERRQ(MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY)); 87 CHKERRQ(MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY)); 88 fclose(Afile); 89 } 90 91 CHKERRQ(PetscOptionsGetString(NULL,NULL,"-rhs",rhs,sizeof(rhs),&flg_b)); 92 if (flg_b) { 93 CHKERRQ(VecCreate(PETSC_COMM_SELF,&b)); 94 CHKERRQ(VecSetSizes(b,PETSC_DECIDE,n)); 95 CHKERRQ(VecSetFromOptions(b)); 96 CHKERRQ(PetscPrintf(PETSC_COMM_SELF,"\n Read rhs in ascii format ...\n")); 97 CHKERRQ(PetscFOpen(PETSC_COMM_SELF,rhs,"r",&bfile)); 98 for (i=0; i<n; i++) { 99 PetscCheckFalse(fscanf(bfile,"%d %le\n",&dummy,(double*)&val) != 2,PETSC_COMM_SELF,PETSC_ERR_FILE_UNEXPECTED,"Badly formatted input file"); 100 CHKERRQ(VecSetValues(b,1,&i,&val,INSERT_VALUES)); 101 } 102 CHKERRQ(VecAssemblyBegin(b)); 103 CHKERRQ(VecAssemblyEnd(b)); 104 fclose(bfile); 105 } 106 107 CHKERRQ(PetscOptionsGetString(NULL,NULL,"-solu",solu,sizeof(solu),&flg_u)); 108 if (flg_u) { 109 CHKERRQ(VecCreate(PETSC_COMM_SELF,&u)); 110 CHKERRQ(VecSetSizes(u,PETSC_DECIDE,n)); 111 CHKERRQ(VecSetFromOptions(u)); 112 CHKERRQ(PetscPrintf(PETSC_COMM_SELF,"\n Read exact solution in ascii format ...\n")); 113 CHKERRQ(PetscFOpen(PETSC_COMM_SELF,solu,"r",&ufile)); 114 for (i=0; i<n; i++) { 115 PetscCheckFalse(fscanf(ufile,"%d %le\n",&dummy,(double*)&val) != 2,PETSC_COMM_SELF,PETSC_ERR_FILE_UNEXPECTED,"Badly formatted input file"); 116 CHKERRQ(VecSetValues(u,1,&i,&val,INSERT_VALUES)); 117 } 118 CHKERRQ(VecAssemblyBegin(u)); 119 CHKERRQ(VecAssemblyEnd(u)); 120 fclose(ufile); 121 } 122 123 /* Write matrix, rhs and exact solution in Petsc binary file */ 124 CHKERRQ(PetscPrintf(PETSC_COMM_SELF,"\n Write matrix in binary to 'matrix.dat' ...\n")); 125 CHKERRQ(PetscViewerBinaryOpen(PETSC_COMM_SELF,"matrix.dat",FILE_MODE_WRITE,&view)); 126 CHKERRQ(MatView(A,view)); 127 128 if (flg_b) { /* Write rhs in Petsc binary file */ 129 CHKERRQ(PetscPrintf(PETSC_COMM_SELF,"\n Write rhs in binary to 'matrix.dat' ...\n")); 130 CHKERRQ(VecView(b,view)); 131 } 132 if (flg_u) { 133 CHKERRQ(PetscPrintf(PETSC_COMM_SELF,"\n Write exact solution in binary to 'matrix.dat' ...\n")); 134 CHKERRQ(VecView(u,view)); 135 } 136 137 /* Check accuracy of the data */ 138 if (flg_A & flg_b & flg_u) { 139 CHKERRQ(VecDuplicate(u,&u_tmp)); 140 CHKERRQ(MatMult(A,u,u_tmp)); 141 CHKERRQ(VecAXPY(u_tmp,-1.0,b)); 142 CHKERRQ(VecNorm(u_tmp,NORM_2,&res_norm)); 143 CHKERRQ(PetscPrintf(PETSC_COMM_SELF,"\n Accuracy of the reading data: | b - A*u |_2 : %g \n",res_norm)); 144 CHKERRQ(VecDestroy(&u_tmp)); 145 } 146 147 CHKERRQ(MatDestroy(&A)); 148 if (flg_b) CHKERRQ(VecDestroy(&b)); 149 if (flg_u) CHKERRQ(VecDestroy(&u)); 150 CHKERRQ(PetscViewerDestroy(&view)); 151 CHKERRQ(PetscFinalize()); 152 return 0; 153 } 154 155 /*TEST 156 157 build: 158 requires: !defined(PETSC_USE_64BIT_INDICES) double !complex 159 160 test: 161 requires: datafilespath 162 args: -Ain ${DATAFILESPATH}/matrices/indefinite/afiro_A.dat 163 164 TEST*/ 165