2: #include <petscmat.h>

  4: static char help[] = "Read in a Symmetric matrix in MatrixMarket format (only the lower triangle). \n\
  5:   Assemble it to a PETSc sparse SBAIJ (upper triangle) matrix. \n\
  6:   Write it in a AIJ matrix (entire matrix) to a file. \n\
  7:   Input parameters are:            \n\
  8:     -fin <filename> : input file   \n\
  9:     -fout <filename> : output file \n\n";

 13: int main(int argc,char **args)
 14: {
 15:   Mat            A;
 16:   char           filein[PETSC_MAX_PATH_LEN],fileout[PETSC_MAX_PATH_LEN],buf[PETSC_MAX_PATH_LEN];
 17:   PetscInt       i,m,n,nnz;
 19:   PetscMPIInt    size;
 20:   PetscScalar    *VAL,zero=0.0;
 21:   FILE           *file;
 22:   PetscViewer    view;
 23:   int            *I,*J,*rownz;

 25:   PetscInitialize(&argc,&args,(char*)0,help);
 26: #if defined(PETSC_USE_COMPLEX)
 27:   SETERRQ(PETSC_COMM_WORLD,1,"This example does not work with complex numbers");
 28: #else
 29:   MPI_Comm_size(PETSC_COMM_WORLD,&size);
 30:   if (size > 1) SETERRQ(PETSC_COMM_WORLD,1,"Uniprocessor Example only\n");

 32:   /* Read in matrix and RHS */
 33:   PetscOptionsGetString(NULL,"-fin",filein,PETSC_MAX_PATH_LEN,NULL);
 34:   PetscFOpen(PETSC_COMM_SELF,filein,"r",&file);

 36:   /* process header with comments */
 37:   do fgets(buf,PETSC_MAX_PATH_LEN-1,file);
 38:   while (buf[0] == '%');

 40:   /* The first non-comment line has the matrix dimensions */
 41:   sscanf(buf,"%d %d %d\n",&m,&n,&nnz);
 42:   PetscPrintf (PETSC_COMM_SELF,"m = %d, n = %d, nnz = %d\n",m,n,nnz);

 44:   /* reseve memory for matrices */
 45:   PetscMalloc4(nnz,&I,nnz,&J,nnz,&VAL,m,&rownz);
 46:   for (i=0; i<m; i++) rownz[i] = 1; /* add 0.0 to diagonal entries */

 48:   for (i=0; i<nnz; i++) {
 49:     fscanf(file,"%d %d %le\n",&I[i],&J[i],(double*)&VAL[i]);
 50:     if (ierr == EOF) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_USER,"i=%d, reach EOF\n",i);
 51:     I[i]--; J[i]--;    /* adjust from 1-based to 0-based */
 52:     rownz[J[i]]++;
 53:   }
 54:   fclose(file);
 55:   PetscPrintf(PETSC_COMM_SELF,"Read file completes.\n");

 57:   /* Creat and asseble SBAIJ matrix */
 58:   MatCreate(PETSC_COMM_SELF,&A);
 59:   MatSetType(A,MATSBAIJ);
 60:   MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,m,n);
 61:   MatSetFromOptions(A);
 62:   MatSeqSBAIJSetPreallocation(A,1,0,rownz);
 63:   MatSetUp(A);

 65:   /* Add zero to diagonals, in case the matrix missing diagonals */
 66:   for (i=0; i<m; i++){
 67:     MatSetValues(A,1,&i,1,&i,&zero,INSERT_VALUES);
 68:   }
 69:   for (i=0; i<nnz; i++) {
 70:     MatSetValues(A,1,&J[i],1,&I[i],&VAL[i],INSERT_VALUES);
 71:   }
 72:   MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);
 73:   MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);
 74:   PetscPrintf(PETSC_COMM_SELF,"Assemble SBAIJ matrix completes.\n");

 76:   /* Write out matrix in AIJ format */
 77:   PetscOptionsGetString(NULL,"-fout",fileout,PETSC_MAX_PATH_LEN,NULL);
 78:   PetscViewerBinaryOpen(PETSC_COMM_WORLD,fileout,FILE_MODE_WRITE,&view);
 79:   MatView(A,view);
 80:   PetscViewerDestroy(&view);

 82:   PetscFree4(I,J,VAL,rownz);
 83:   MatDestroy(&A);
 84:   PetscFinalize();
 85: #endif
 86:   return 0;
 87: }