Actual source code: ex78.c
petsc-3.12.4 2020-02-04
Report Typos and Errors

      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";

     10: /*
     11: Example: ./ex78 -Ain Ain -rhs rhs -solu solu -noshift -mat_view
     12:  with the datafiles in the followig 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

     23: rhs
     24: ---
     25: 0 3.0
     26: 1 12.0
     27: 2 6.0

     29: solu
     30: ----
     31: 0 1.0
     32: 0 1.0
     33: 0 1.0
     34: */

     37:  #include <petscmat.h>

     39: int main(int argc,char **args)
     40: {
     41:   Mat            A = NULL;
     42:   Vec            b,u = NULL,u_tmp;
     43:   char           Ain[PETSC_MAX_PATH_LEN],rhs[PETSC_MAX_PATH_LEN],solu[PETSC_MAX_PATH_LEN];
     45:   int            m,n = 0,nz,dummy; /* these are fscaned so kept as int */
     46:   PetscInt       i,col,row,shift = 1,sizes[3],nsizes;
     47:   PetscScalar    val;
     48:   PetscReal      res_norm;
     49:   FILE           *Afile,*bfile,*ufile;
     50:   PetscViewer    view;
     51:   PetscBool      flg_A,flg_b,flg_u,flg;
     52:   PetscMPIInt    size;

     54:   PetscInitialize(&argc,&args,(char*)0,help);if (ierr) return ierr;
     55:   MPI_Comm_size(PETSC_COMM_WORLD,&size);
     56:   if (size != 1) SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_SUP,"This is a uniprocessor example only!");

     58:   /* Read in matrix, rhs and exact solution from ascii files */
     59:   PetscOptionsGetString(NULL,NULL,"-Ain",Ain,PETSC_MAX_PATH_LEN,&flg_A);
     60:   PetscOptionsHasName(NULL,NULL,"-noshift",&flg);
     61:   if (flg) shift = 0;
     62:   if (flg_A) {
     63:     PetscPrintf(PETSC_COMM_SELF,"\n Read matrix in ascii format ...\n");
     64:     PetscFOpen(PETSC_COMM_SELF,Ain,"r",&Afile);
     65:     nsizes = 3;
     66:     PetscOptionsGetIntArray(NULL,NULL,"-nosizesinfile",sizes,&nsizes,&flg);
     67:     if (flg) {
     68:       if (nsizes != 3) SETERRQ(PETSC_COMM_WORLD,1,"Must pass in three m,n,nz as arguments for -nosizesinfile");
     69:       m  = sizes[0];
     70:       n  = sizes[1];
     71:       nz = sizes[2];
     72:     } else {
     73:       if (fscanf(Afile,"%d %d %d\n",&m,&n,&nz) != 3)  SETERRQ(PETSC_COMM_SELF,1,"Badly formatted input file\n");
     74:     }
     75:     PetscPrintf(PETSC_COMM_SELF,"m: %d, n: %d, nz: %d \n", m,n,nz);
     76:     if (m != n) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_SIZ, "Number of rows, cols must be same for this example\n");
     77:     MatCreate(PETSC_COMM_SELF,&A);
     78:     MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,m,n);
     79:     MatSetFromOptions(A);
     80:     MatSeqAIJSetPreallocation(A,nz/m,NULL);
     81:     MatSetOption(A,MAT_NEW_NONZERO_ALLOCATION_ERR,PETSC_FALSE);

     83:     for (i=0; i<nz; i++) {
     84:       if (fscanf(Afile,"%d %d %le\n",&row,&col,(double*)&val) != 3) SETERRQ(PETSC_COMM_SELF,1,"Badly formatted input file\n");
     85:       row -= shift; col -= shift;  /* set index set starts at 0 */
     86:       MatSetValues(A,1,&row,1,&col,&val,INSERT_VALUES);
     87:     }
     88:     MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);
     89:     MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);
     90:     fclose(Afile);
     91:   }

     93:   PetscOptionsGetString(NULL,NULL,"-rhs",rhs,PETSC_MAX_PATH_LEN,&flg_b);
     94:   if (flg_b) {
     95:     VecCreate(PETSC_COMM_SELF,&b);
     96:     VecSetSizes(b,PETSC_DECIDE,n);
     97:     VecSetFromOptions(b);
     98:     PetscPrintf(PETSC_COMM_SELF,"\n Read rhs in ascii format ...\n");
     99:     PetscFOpen(PETSC_COMM_SELF,rhs,"r",&bfile);
    100:     for (i=0; i<n; i++) {
    101:       if (fscanf(bfile,"%d %le\n",&dummy,(double*)&val) != 2) SETERRQ(PETSC_COMM_SELF,1,"Badly formatted input file\n");
    102:       VecSetValues(b,1,&i,&val,INSERT_VALUES);
    103:     }
    104:     VecAssemblyBegin(b);
    105:     VecAssemblyEnd(b);
    106:     fclose(bfile);
    107:   }

    109:   PetscOptionsGetString(NULL,NULL,"-solu",solu,PETSC_MAX_PATH_LEN,&flg_u);
    110:   if (flg_u) {
    111:     VecCreate(PETSC_COMM_SELF,&u);
    112:     VecSetSizes(u,PETSC_DECIDE,n);
    113:     VecSetFromOptions(u);
    114:     PetscPrintf(PETSC_COMM_SELF,"\n Read exact solution in ascii format ...\n");
    115:     PetscFOpen(PETSC_COMM_SELF,solu,"r",&ufile);
    116:     for (i=0; i<n; i++) {
    117:       if (fscanf(ufile,"%d  %le\n",&dummy,(double*)&val) != 2) SETERRQ(PETSC_COMM_SELF,1,"Badly formatted input file\n");
    118:       VecSetValues(u,1,&i,&val,INSERT_VALUES);
    119:     }
    120:     VecAssemblyBegin(u);
    121:     VecAssemblyEnd(u);
    122:     fclose(ufile);
    123:   }

    125:   /* Write matrix, rhs and exact solution in Petsc binary file */
    126:   PetscPrintf(PETSC_COMM_SELF,"\n Write matrix in binary to 'matrix.dat' ...\n");
    127:   PetscViewerBinaryOpen(PETSC_COMM_SELF,"matrix.dat",FILE_MODE_WRITE,&view);
    128:   MatView(A,view);

    130:   if (flg_b) { /* Write rhs in Petsc binary file */
    131:     PetscPrintf(PETSC_COMM_SELF,"\n Write rhs in binary to 'matrix.dat' ...\n");
    132:     VecView(b,view);
    133:   }
    134:   if (flg_u) {
    135:     PetscPrintf(PETSC_COMM_SELF,"\n Write exact solution in binary to 'matrix.dat' ...\n");
    136:     VecView(u,view);
    137:   }

    139:   /* Check accuracy of the data */
    140:   if (flg_A & flg_b & flg_u) {
    141:     VecDuplicate(u,&u_tmp);
    142:     MatMult(A,u,u_tmp);
    143:     VecAXPY(u_tmp,-1.0,b);
    144:     VecNorm(u_tmp,NORM_2,&res_norm);
    145:     PetscPrintf(PETSC_COMM_SELF,"\n Accuracy of the reading data: | b - A*u |_2 : %g \n",res_norm);
    146:     VecDestroy(&u_tmp);
    147:   }

    149:   MatDestroy(&A);
    150:   if (flg_b) {VecDestroy(&b);}
    151:   if (flg_u) {VecDestroy(&u);}
    152:   PetscViewerDestroy(&view);
    153:   PetscFinalize();
    154:   return ierr;
    155: }

    157: /*TEST

    159:    build:
    160:       requires:  !define(PETSC_USE_64BIT_INDICES) double !complex datafilespath

    162:    test:
    163:       args: -Ain ${DATAFILESPATH}/matrices/indefinite/afiro_A.dat

    165: TEST*/