2: static char help[] = "Reads a matrix from PETSc binary file. Use for view or investigating matrix data structure. \n\n";
  3: /*
  4:  Example:
  5:       ./ex16 -f <matrix file> -a_mat_view draw -draw_pause -1
  6:       ./ex16 -f <matrix file> -a_mat_view ascii::ascii_info
  7:  */

  9: #include <petscmat.h>
 10: int main(int argc, char **args)
 11: {
 12:   Mat                A, Asp;
 13:   PetscViewer        fd;                       /* viewer */
 14:   char               file[PETSC_MAX_PATH_LEN]; /* input file name */
 15:   PetscInt           m, n, rstart, rend;
 16:   PetscBool          flg;
 17:   PetscInt           row, ncols, j, nrows, nnzA = 0, nnzAsp = 0;
 18:   const PetscInt    *cols;
 19:   const PetscScalar *vals;
 20:   PetscReal          norm, percent, val, dtol = 1.e-16;
 21:   PetscMPIInt        rank;
 22:   MatInfo            matinfo;
 23:   PetscInt           Dnnz, Onnz;

 26:   PetscInitialize(&argc, &args, (char *)0, help);
 27:   MPI_Comm_rank(PETSC_COMM_WORLD, &rank);

 29:   /* Determine files from which we read the linear systems. */
 30:   PetscOptionsGetString(NULL, NULL, "-f", file, sizeof(file), &flg);

 33:   /* Open binary file.  Note that we use FILE_MODE_READ to indicate
 34:      reading from this file. */
 35:   PetscViewerBinaryOpen(PETSC_COMM_WORLD, file, FILE_MODE_READ, &fd);

 37:   /* Load the matrix; then destroy the viewer. */
 38:   MatCreate(PETSC_COMM_WORLD, &A);
 39:   MatSetOptionsPrefix(A, "a_");
 40:   MatSetFromOptions(A);
 41:   MatLoad(A, fd);
 42:   PetscViewerDestroy(&fd);
 43:   MatGetSize(A, &m, &n);
 44:   MatGetInfo(A, MAT_LOCAL, &matinfo);
 45:   /*printf("matinfo.nz_used %g\n",matinfo.nz_used);*/

 47:   /* Get a sparse matrix Asp by dumping zero entries of A */
 48:   MatCreate(PETSC_COMM_WORLD, &Asp);
 49:   MatSetSizes(Asp, m, n, PETSC_DECIDE, PETSC_DECIDE);
 50:   MatSetOptionsPrefix(Asp, "asp_");
 51:   MatSetFromOptions(Asp);
 52:   Dnnz = (PetscInt)matinfo.nz_used / m + 1;
 53:   Onnz = Dnnz / 2;
 54:   printf("Dnnz %d %d\n", Dnnz, Onnz);
 55:   MatSeqAIJSetPreallocation(Asp, Dnnz, NULL);
 56:   MatMPIAIJSetPreallocation(Asp, Dnnz, NULL, Onnz, NULL);
 57:   /* The allocation above is approximate so we must set this option to be permissive.
 58:    * Real code should preallocate exactly. */
 59:   MatSetOption(Asp, MAT_NEW_NONZERO_LOCATION_ERR, PETSC_FALSE);

 61:   /* Check zero rows */
 62:   MatGetOwnershipRange(A, &rstart, &rend);
 63:   nrows = 0;
 64:   for (row = rstart; row < rend; row++) {
 65:     MatGetRow(A, row, &ncols, &cols, &vals);
 66:     nnzA += ncols;
 67:     norm = 0.0;
 68:     for (j = 0; j < ncols; j++) {
 69:       val = PetscAbsScalar(vals[j]);
 70:       if (norm < val) norm = norm;
 71:       if (val > dtol) {
 72:         MatSetValues(Asp, 1, &row, 1, &cols[j], &vals[j], INSERT_VALUES);
 73:         nnzAsp++;
 74:       }
 75:     }
 76:     if (!norm) nrows++;
 77:     MatRestoreRow(A, row, &ncols, &cols, &vals);
 78:   }
 79:   MatAssemblyBegin(Asp, MAT_FINAL_ASSEMBLY);
 80:   MatAssemblyEnd(Asp, MAT_FINAL_ASSEMBLY);

 82:   percent = (PetscReal)nnzA * 100 / (m * n);
 83:   PetscPrintf(PETSC_COMM_SELF, " [%d] Matrix A local size %d,%d; nnzA %d, %g percent; No. of zero rows: %d\n", rank, m, n, nnzA, percent, nrows);
 84:   percent = (PetscReal)nnzAsp * 100 / (m * n);
 85:   PetscPrintf(PETSC_COMM_SELF, " [%d] Matrix Asp nnzAsp %d, %g percent\n", rank, nnzAsp, percent);

 87:   /* investigate matcoloring for Asp */
 88:   PetscBool Asp_coloring = PETSC_FALSE;
 89:   PetscOptionsHasName(NULL, NULL, "-Asp_color", &Asp_coloring);
 90:   if (Asp_coloring) {
 91:     MatColoring   mc;
 92:     ISColoring    iscoloring;
 93:     MatFDColoring matfdcoloring;
 94:     PetscPrintf(PETSC_COMM_WORLD, " Create coloring of Asp...\n");
 95:     MatColoringCreate(Asp, &mc);
 96:     MatColoringSetType(mc, MATCOLORINGSL);
 97:     MatColoringSetFromOptions(mc);
 98:     MatColoringApply(mc, &iscoloring);
 99:     MatColoringDestroy(&mc);
100:     MatFDColoringCreate(Asp, iscoloring, &matfdcoloring);
101:     MatFDColoringSetFromOptions(matfdcoloring);
102:     MatFDColoringSetUp(Asp, iscoloring, matfdcoloring);
103:     /*MatFDColoringView(matfdcoloring,PETSC_VIEWER_STDOUT_WORLD);*/
104:     ISColoringDestroy(&iscoloring);
105:     MatFDColoringDestroy(&matfdcoloring);
106:   }

108:   /* Write Asp in binary for study - see ~petsc/src/mat/tests/ex124.c */
109:   PetscBool Asp_write = PETSC_FALSE;
110:   PetscOptionsHasName(NULL, NULL, "-Asp_write", &Asp_write);
111:   if (Asp_write) {
112:     PetscViewer viewer;
113:     PetscPrintf(PETSC_COMM_SELF, "Write Asp into file Asp.dat ...\n");
114:     PetscViewerBinaryOpen(PETSC_COMM_WORLD, "Asp.dat", FILE_MODE_WRITE, &viewer);
115:     MatView(Asp, viewer);
116:     PetscViewerDestroy(&viewer);
117:   }

119:   MatDestroy(&A);
120:   MatDestroy(&Asp);
121:   PetscFinalize();
122:   return 0;
123: }